Skip to content

Commit

Permalink
Merge remote-tracking branch 'origin/master'
Browse files Browse the repository at this point in the history 
# Conflicts:
#	src/hic/tools/utils/norm/RU/RUCode.java
  • Loading branch information
@sa501428
sa501428 committed Sep 1, 2022
2 parents e1af332 + 40fa453 commit 3504ae8
Showing 1 changed file with 70 additions and 44 deletions.
114 changes: 70 additions & 44 deletions src/hic/tools/utils/norm/RU/RUCode.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -24,12 +24,13 @@

package hic.tools.utils.norm.RU;

public class RUCode {
import java.util.Arrays;

public class RUCode() {
public static int balance(long m, double[] b, double[] report, double[] allIters, float tol, float pppp,
int maxiter, double del, double dp, int totIter, int threads, int k) {

/*
int p;
int n0;
double low;
int[] bad;
Expand All @@ -49,53 +50,63 @@ public static int balance(long m, double[] b, double[] report, double[] allIters

double[] space = new double[threads];

for (int p = 0; p < threads; p++) {
space[p] = new double[k];
for (int c = 0; c < threads; c++) {
space[c] = new double[k];
}

double perc = pppp;

for (int p = 0; p < k; p++){
one[p] = 1.0;
for (int g = 0; g < k; g++){
one[g] = 1.0;
}
for (int p = 0; p < k; p++){
bad[p] = 0;

for (int g = 0; g < k; g++){
bad[g] = 0;
}


double[] nz = new double[k];
for (int p = 0; p < k; p++){
nz[p] = 0;

for (int g = 0; g < k; g++){
nz[g] = 0;
}

for (int p = 0; p < m; p++) {
for (p = 0; p < m; p++) {
nz[i[p]] += 1.0;
if (i[p] != j[p]) nz[j[p]] += 1.0;
}

n0 = 0;
for (p = 0; p < k; p++) if (nz[p] > 0) r0[n0++] = nz[p];
for (int g = 0; g < k; g++) if (nz[g] > 0) r0[n0++] = nz[g];
Arrays.sort(r0);

lind = (int) (n0 * perc + 0.5);
if (lind < 0) lind = 0;
low = r0[lind];
// find rows which are identically 0 and the perc proportion of rows with the lowest row sums and exclude them
for (p = 0; p < k; p++) if (nz[p] < low) bad[p] = 1;
for (int g = 0; g < k; g++) if (nz[g] < low) bad[g] = 1;
// if bad[p] is 1 we remove row p
utmvMul(i, j, x, m, one, k, row0, threads, space);
for (p = 0; p < k; p++) row[p] = row0[p];
for (int g = 0; g < k; g++){
row[g] = row0[g];
}

for (p = 0; p < k; p++) one[p] = 1.0 - bad[p];
for (p = 0; p < k; p++) if (bad[p] == 1) row[p] = 1.0;
for (p = 0; p < k; p++) b[p] = one[p] / Math.sqrt(row[p]);
for (int g = 0; g < k; g++){
one[g] = 1.0 - bad[g];
}
for (int g = 0; g < k; g++) if (bad[g] == 1) row[g] = 1.0;
for (int g = 0; g < k; g++){
b[g] = one[g] / Math.sqrt(row[g]);
}

// start iterations
// row is the current rows sum; dr and dc are the current rows and columns scaling vectors
double ber = 10.0 * (1.0 + tol);
double err = 10.0 * (1.0 + tol);
int iter = 0;
for (p = 0; p < k; p++) current[p] = b[p];
for (int g = 0; g < k; g++){
current[g] = b[g];
}
int fail; // checks whether the convergence rate seems to be good enough; 0 if yes, 1 if no
int nerr = 0;

Expand All @@ -108,15 +119,21 @@ public static int balance(long m, double[] b, double[] report, double[] allIters
fail = 1;

utmvMul(i, j, x, m, b, k, row, threads, space);
for (p = 0; p < k; p++) row[p] *= b[p];
for (p = 0; p < k; p++) if (bad[p] == 1) row[p] = 1.0;
for (p = 0; p < k; p++) b[p] *= one[p] / Math.sqrt(row[p]);
for (int g = 0; g < k; g++){
row[g] *= b[g];
}
for (int g = 0; g < k; g++) if (bad[g] == 1){
row[g] = 1.0;
}
for (int g = 0; g < k; g++){
b[g] *= one[g] / Math.sqrt(row[g]);
}
// calculate the current error
ber = 0;
for (p = 0; p < k; p++) {
if (bad[p] == 1) continue;
if (Math.abs((b[p] - current[p]) / (b[p] + current[p])) > ber)
ber = Math.abs((b[p] - current[p]) / (b[p] + current[p]));
for (int g = 0; g < k; g++) {
if (bad[g] == 1) continue;
if (Math.abs((b[g] - current[g]) / (b[g] + current[g])) > ber)
ber = Math.abs((b[g] - current[g]) / (b[g] + current[g]));
}
report[all_iters - 1] = ber;
allIters[all_iters - 1] = iter;
Expand All @@ -126,20 +143,22 @@ public static int balance(long m, double[] b, double[] report, double[] allIters
if (iter % 10 == 0) {
utmvMul(i, j, x, m, b, k, row, threads, space);
err = 0;
for (p = 0; p < k; p++) {
if (bad[p] == 1) continue;
if (err < Math.abs(row[p] * b[p] - one[p]))
err = Math.abs(row[p] * b[p] - one[p]);
for (int g = 0; g < k; g++) {
if (bad[g] == 1) continue;
if (err < Math.abs(row[g] * b[g] - one[g]))
err = Math.abs(row[g] * b[g] - one[g]);
}
errors[nerr++] = err;
}
for (p = 0; p < k; p++) current[p] = b[p];
for (int g = 0; g < k; g++){
current[g] = b[g];
}
// check whether convergence rate is satisfactory
// if less than 5 iterations (so less than 5 errors) and less than 2 row sums errors, there is nothing to chek
if ((ber < tol) && (nerr < 2 || (nerr >= 2 && errors[nerr - 1] < 0.5 * errors[nerr - 2]))) continue;
// otherwise check
if (iter > 5) {
for (p = 1; p <= 5; p++) if (report[all_iters - p] * (1.0 + del) < report[all_iters - p - 1]) fail = 0;
for (int g = 1; g <= 5; g++) if (report[all_iters - g] * (1.0 + del) < report[all_iters - g - 1]) fail = 0;
if (nerr >= 2 && errors[nerr - 1] > 0.75 * errors[nerr - 2]) fail = 1;
if (iter >= maxiter) fail = 1;

Expand All @@ -150,10 +169,10 @@ public static int balance(long m, double[] b, double[] report, double[] allIters
nerr = 0;
lind = (int) (n0 * perc + 0.5);
low = r0[lind];
for (p = 0; p < k; p++) {
if (nz[p] < low) {
bad[p] = 1;
one[p] = 0;
for (int g = 0; g < k; g++) {
if (nz[g] < low) {
bad[g] = 1;
one[g] = 0;
}
}

Expand All @@ -162,11 +181,21 @@ public static int balance(long m, double[] b, double[] report, double[] allIters
// if the current error is larger than 5 iteration ago start from scratch, otherwise continue from the current
// position
if (report[all_iters - 1] > report[all_iters - 6]) {
for (p = 0; p < k; p++) one[p] = 1.0 - bad[p];
for (p = 0; p < k; p++) row0[p] *= one[p];
for (p = 0; p < k; p++) row[p] = row0[p];
for (p = 0; p < k; p++) current[p] = b[p];
} else for (p = 0; p < k; p++) b[p] *= (1.0 - bad[p]);
for (p = 0; p < k; p++){
one[p] = 1.0 - bad[p];
}
for (p = 0; p < k; p++){
row0[p] *= one[p];
}
for (p = 0; p < k; p++){
row[p] = row0[p];
}
for (p = 0; p < k; p++){
current[p] = b[p];
}
} else for (p = 0; p < k; p++){
b[p] *= (1.0 - bad[p]);
}
iter = 0;

}
Expand All @@ -190,9 +219,6 @@ public static int balance(long m, double[] b, double[] report, double[] allIters
for (p = 0; p < k; p++) if (bad[p] == 1) b[p] = Float.NaN;

return (iter);
*/
return 0;
}
}

0 comments on commit 3504ae8

Please sign in to comment.