Skip to content
Snippets Groups Projects
Commit dc7ff81a authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

The sorted interactions test case can now be run for each of the three possible orientations.

parent 785f596e
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
examples/plot_sorted.py
+ 111
− 107
View file @ dc7ff81a
......@@ -36,113 +36,117 @@ import os
print "Plotting..."
# Read Quickshed accelerations
data=loadtxt("interaction_dump.dat")
id = data[:,0]
posx = data[:,2]
accx_u=data[:,5]
accy_u=data[:,6]
accz_u=data[:,7]
accx_s=data[:,8]
accy_s=data[:,9]
accz_s=data[:,10]
# Build error ------------------------------------------------
errx_s = (accx_s - accx_u )/abs(accx_u)
erry_s = (accy_s - accy_u )/abs(accy_u)
errz_s = (accz_s - accz_u )/abs(accz_u)
e_errx_s = errx_s#[abs(errx_s) > 0.001]
e_erry_s = erry_s#[abs(erry_s) > 0.001]
e_errz_s = errz_s#[abs(errz_s) > 0.001]
# Statistics
meanx_s = mean(errx_s[abs(errx_s) < 0.1])
stdx_s = std(errx_s[abs(errx_s) < 0.1])
meany_s = mean(erry_s[abs(erry_s) < 0.1])
stdy_s = std(erry_s[abs(erry_s) < 0.1])
meanz_s = mean(errz_s[abs(errz_s) < 0.1])
stdz_s = std(errz_s[abs(errz_s) < 0.1])
# Plot -------------------------------------------------------
figure(frameon=True)
subplot(311, title="Acceleration along X")
#plot(id[abs(errx_s) > 0.001], e_errx_s , 'rx')
plot(posx, e_errx_s , 'rx')
#text(id[-1], 0.18, "B-H: $%5.3f\\pm%5.3f$\n QuickShed: $%5.3f\\pm%5.3f$"%(meanx_bh, stdx_bh, meanx_new, stdx_new), backgroundcolor="w", va="top", ha="right" )
ylim(-0.2, 0.2)
#xlim(0, size(id)-1)
grid()
subplot(312, title="Acceleration along Y")
#plot(id[abs(erry_s) > 0.001], e_erry_s , 'rx')
plot(posx, e_erry_s , 'rx')
#text(id[-1], 0.18, "B-H: $%5.3f\\pm%5.3f$\n QuickShed: $%5.3f\\pm%5.3f$"%(meany_bh, stdy_bh, meany_new, stdy_new), backgroundcolor="w", va="top", ha="right" )
ylim(-0.2, 0.2)
#xlim(0, size(id)-1)
grid()
subplot(313, title="Acceleration along Z")
#plot(id[abs(errz_s) > 0.001], e_errz_s , 'rx', label="Sorted")
plot(posx, e_errz_s , 'rx', label="Sorted")
#text(id[-1], 0.18, "B-H: $%5.3f\\pm%5.3f$\n QuickShed: $%5.3f\\pm%5.3f$"%(meanz_bh, stdz_bh, meanz_new, stdz_new), backgroundcolor="w", va="top", ha="right" )
legend(loc="upper right")
ylim(-0.2, 0.2)
#xlim(0, size(id)-1)
grid()
savefig("accelerations.png")
# Plot -------------------------------------------------------
bins = linspace(-3, 3, 10000)
e_errx_s = errx_s[(abs(errx_s) > 0.001) & (abs(errx_s) < 1.)]
e_erry_s = erry_s[(abs(erry_s) > 0.001) & (abs(erry_s) < 1.)]
e_errz_s = errz_s[(abs(errz_s) > 0.001) & (abs(errz_s) < 1.)]
figure(frameon=True)
subplot(311, title="Acceleration along X")#, yscale='log')
hist(e_errx_s, bins=bins, normed=1, histtype='step', rwidth=0.01, color='r', label="Sorted")
legend(loc="upper right")
xlim(-0.12, 0.12)
subplot(312, title="Acceleration along Y")
hist(e_erry_s, bins=bins, normed=1, histtype='step', rwidth=0.01, color='r')
xlim(-0.12, 0.12)
subplot(313, title="Acceleration along Z")
hist(e_errz_s, bins=bins, normed=1, histtype='step', rwidth=0.01, color='r')
xlim(-0.12, 0.12)
savefig("histogram.png")
print "E_error for sorted: ( x= %4.3e"%std(e_errx_s), "y= %4.3e"%std(e_erry_s), "z= %4.3e"%std(e_errz_s), ") Combined YZ = %4.3e"%(( std(e_erry_s) + std(e_errz_s) )/2.)
#print std(e_errx_s), std(e_erry_s), std(e_errz_s)
#print "Min for sorted: ( x= %4.3e"%min(errx_s), "y= %4.3e"%min(erry_s), "z= %4.3e"%min(errz_s), ") Combined YZ = %4.3e"%(min( min(erry_s), min(errz_s) ))
#print "Max for sorted: ( x= %4.3e"%max(errx_s), "y= %4.3e"%max(erry_s), "z= %4.3e"%max(errz_s), ") Combined YZ = %4.3e"%(max( max(erry_s), max(errz_s) ))
names = ["side", "edge", "corner"]
for orientation in range( 3 ):
# Read Quickshed accelerations
data=loadtxt( "interaction_dump_%d.dat"%orientation )
id = data[:,0]
pos = data[:,2]
accx_u=data[:,6]
accy_u=data[:,7]
accz_u=data[:,8]
accx_s=data[:,9]
accy_s=data[:,10]
accz_s=data[:,11]
# Build error ------------------------------------------------
errx_s = (accx_s - accx_u )/abs(accx_u)
erry_s = (accy_s - accy_u )/abs(accy_u)
errz_s = (accz_s - accz_u )/abs(accz_u)
e_errx_s = errx_s#[abs(errx_s) > 0.001]
e_erry_s = erry_s#[abs(erry_s) > 0.001]
e_errz_s = errz_s#[abs(errz_s) > 0.001]
# Statistics
meanx_s = mean(errx_s[abs(errx_s) < 0.1])
stdx_s = std(errx_s[abs(errx_s) < 0.1])
meany_s = mean(erry_s[abs(erry_s) < 0.1])
stdy_s = std(erry_s[abs(erry_s) < 0.1])
meanz_s = mean(errz_s[abs(errz_s) < 0.1])
stdz_s = std(errz_s[abs(errz_s) < 0.1])
# Plot -------------------------------------------------------
figure(frameon=True)
subplot(311, title="Acceleration along X")
#plot(id[abs(errx_s) > 0.001], e_errx_s , 'rx')
plot(pos, e_errx_s , 'rx')
#text(id[-1], 0.18, "B-H: $%5.3f\\pm%5.3f$\n QuickShed: $%5.3f\\pm%5.3f$"%(meanx_bh, stdx_bh, meanx_new, stdx_new), backgroundcolor="w", va="top", ha="right" )
ylim(-0.2, 0.2)
#xlim(0, size(id)-1)
grid()
subplot(312, title="Acceleration along Y")
#plot(id[abs(erry_s) > 0.001], e_erry_s , 'rx')
plot(pos, e_erry_s , 'rx')
#text(id[-1], 0.18, "B-H: $%5.3f\\pm%5.3f$\n QuickShed: $%5.3f\\pm%5.3f$"%(meany_bh, stdy_bh, meany_new, stdy_new), backgroundcolor="w", va="top", ha="right" )
ylim(-0.2, 0.2)
#xlim(0, size(id)-1)
grid()
subplot(313, title="Acceleration along Z")
#plot(id[abs(errz_s) > 0.001], e_errz_s , 'rx', label="Sorted")
plot(pos, e_errz_s , 'rx', label="Sorted")
#text(id[-1], 0.18, "B-H: $%5.3f\\pm%5.3f$\n QuickShed: $%5.3f\\pm%5.3f$"%(meanz_bh, stdz_bh, meanz_new, stdz_new), backgroundcolor="w", va="top", ha="right" )
legend(loc="upper right")
ylim(-0.2, 0.2)
#xlim(0, size(id)-1)
grid()
savefig("accelerations_%s.png"% names[orientation])
# Plot -------------------------------------------------------
bins = linspace(-3, 3, 10000)
e_errx_s = errx_s[(abs(errx_s) > 0.001) & (abs(errx_s) < 1.)]
e_erry_s = erry_s[(abs(erry_s) > 0.001) & (abs(erry_s) < 1.)]
e_errz_s = errz_s[(abs(errz_s) > 0.001) & (abs(errz_s) < 1.)]
figure(frameon=True)
subplot(311, title="Acceleration along X")#, yscale='log')
hist(e_errx_s, bins=bins, normed=1, histtype='step', rwidth=0.01, color='r', label="Sorted")
legend(loc="upper right")
xlim(-0.12, 0.12)
subplot(312, title="Acceleration along Y")
hist(e_erry_s, bins=bins, normed=1, histtype='step', rwidth=0.01, color='r')
xlim(-0.12, 0.12)
subplot(313, title="Acceleration along Z")
hist(e_errz_s, bins=bins, normed=1, histtype='step', rwidth=0.01, color='r')
xlim(-0.12, 0.12)
savefig("histogram_%s.png"% names[orientation])
print "E_error for sorted: ( x= %4.3e"%std(e_errx_s), "y= %4.3e"%std(e_erry_s), "z= %4.3e"%std(e_errz_s), ") Combined YZ = %4.3e"%(( std(e_erry_s) + std(e_errz_s) )/2.)
#print std(e_errx_s), std(e_erry_s), std(e_errz_s)
#print "Min for sorted: ( x= %4.3e"%min(errx_s), "y= %4.3e"%min(erry_s), "z= %4.3e"%min(errz_s), ") Combined YZ = %4.3e"%(min( min(erry_s), min(errz_s) ))
#print "Max for sorted: ( x= %4.3e"%max(errx_s), "y= %4.3e"%max(erry_s), "z= %4.3e"%max(errz_s), ") Combined YZ = %4.3e"%(max( max(erry_s), max(errz_s) ))
examples/test_bh_sorted.c
+ 59
− 17
View file @ dc7ff81a
......@@ -44,7 +44,7 @@
#define dist_cutoff_ratio 1.5
#define iact_pair_direct iact_pair_direct_sorted
#define ICHECK 2
#define ICHECK -1
#define NO_SANITY_CHECKS
#define NO_COM_AS_TASK
#define COUNTERS
......@@ -2018,8 +2018,10 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) {
/**
* @brief Creates two neighbouring cells wiht N_parts per celland makes them interact using both the
* sorted and unsortde interactions. Outputs then the tow sets of accelerations for accuracy tests.
* @param N_parts Number of particles in each cell
* @param orientation Orientation of the cells ( 0 == side, 1 == edge, 2 == corner )
*/
void test_direct_neighbour( int N_parts ) {
void test_direct_neighbour( int N_parts , int orientation ) {
int k;
struct part *parts;
......@@ -2032,19 +2034,34 @@ void test_direct_neighbour( int N_parts ) {
/* Create random set of particles in both cells */
for (k = 0; k < 2*N_parts; k++) {
parts[k].id = k;
parts[k].x[0] = ((double)rand()) / RAND_MAX + (k >= N_parts ? 1 : 0. );
parts[k].x[0] = ((double)rand()) / RAND_MAX;
if ( k >= N_parts )
parts[k].x[0] += 1.;
parts[k].x[1] = ((double)rand()) / RAND_MAX;
if ( orientation > 0 && k >= N_parts )
parts[k].x[1] += 1.;
parts[k].x[2] = ((double)rand()) / RAND_MAX;
parts[k].mass = 1.;//((double)rand()) / RAND_MAX;
if ( orientation > 1 && k >= N_parts )
parts[k].x[2] += 1.;
parts[k].mass = ((double)rand()) / RAND_MAX;
parts[k].a[0] = 0.0;
parts[k].a[1] = 0.0;
parts[k].a[2] = 0.0;
}
/* Get the cell geometry right */
left.loc[0] = 0.; left.loc[1] = 0.; left.loc[2] = 0.;
left.loc[0] = 0.;
left.loc[1] = 0.;
left.loc[2] = 0.;
left.h = 1.;
right.loc[0] = 1.; right.loc[1] = 0.; right.loc[2] = 0.;
right.loc[0] = 1.;
right.loc[1] = ( orientation > 0 ? 1 : 0 );
right.loc[2] = ( orientation > 1 ? 1 : 0 );
right.h = 1.;
/* Put the particles in the cell */
......@@ -2107,27 +2124,50 @@ void test_direct_neighbour( int N_parts ) {
}
}
/* Sort the particles along axis to simply interpretation */
int compParts(const void* c1, const void* c2) {
if ( ((struct part*)c1)->x[0] < ((struct part*)c2)->x[0] ) return -1;
else if ( ((struct part*)c1)->x[0] == ((struct part*)c2)->x[0] ) return 0;
else return 1;
/* Position along the axis */
double* position;
if ((position = (double *)malloc(sizeof(double) * N_parts * 2)) == NULL)
error("Failed to allocate position buffer.");
for ( k = 0; k < 2*N_parts ; ++k ) {
switch( orientation ) {
case 0:
position[k] = parts[k].x[0] - 1.;
break;
case 1:
position[k] = sqrt( ( parts[k].x[0] * parts[k].x[0] ) + ( parts[k].x[1] * parts[k].x[1] ) ) - sqrt(2.);
break;
case 2:
position[k] = sqrt( ( parts[k].x[0] * parts[k].x[0] ) + ( parts[k].x[1] * parts[k].x[1] ) + ( parts[k].x[2] * parts[k].x[2] ) ) - sqrt(3.);
break;
default:
error("Wrong switch statement");
break;
}
}
//qsort( parts, 2*N_parts, sizeof(struct part), compParts);
/* Now, output everything */
message( "Writing file 'interaction_dump.dat'" );
FILE* file = fopen("interaction_dump.dat", "w");
fprintf(file, "# ID m x y z a_u.x a_u.y a_u.z a_s.x a_s.y a_s.z\n");
char fileName[100];
sprintf( fileName, "interaction_dump_%d.dat", orientation );
message( "Writing file '%s'", fileName );
FILE* file = fopen( fileName , "w" );
fprintf(file, "# ID m r x y z a_u.x a_u.y a_u.z a_s.x a_s.y a_s.z\n");
for (k = 0; k < 2*N_parts; k++) {
fprintf(file, "%d %e %e %e %e %e %e %e %e %e %e\n", parts[k].id, parts[k].mass,
fprintf(file, "%d %e %e %e %e %e %e %e %e %e %e %e\n", parts[k].id, parts[k].mass, position[k],
parts[k].x[0], parts[k].x[1], parts[k].x[2],
parts[k].a_exact[0], parts[k].a_exact[1], parts[k].a_exact[2],
parts[k].a[0], parts[k].a[1], parts[k].a[2]);
}
fclose( file );
#ifdef COUNTERS
message( "Unsorted intereactions: %d" ,count_direct_unsorted );
message( "Sorted intereactions PP: %d" ,count_direct_sorted_pp );
......@@ -2211,7 +2251,9 @@ int main(int argc, char *argv[]) {
message("Interacting 2 neighbouring cells with %d particles per cell", N_parts);
/* Run the test */
test_direct_neighbour( N_parts );
test_direct_neighbour( N_parts, 0 );
test_direct_neighbour( N_parts, 1 );
test_direct_neighbour( N_parts, 2 );
}
else {
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment