Skip to content
Snippets Groups Projects
Commit 37ae4986 authored by Tobias Baetge's avatar Tobias Baetge Committed by Cristian Lalescu
Browse files

cleaning up for merging

parent 7e740857
No related branches found
No related tags found
No related merge requests found
......@@ -35,25 +35,23 @@ import matplotlib.pyplot as plt
import TurTLE
from TurTLE import DNS
def theory_exponential(t,initial_value,drag_coeff):
return initial_value*np.exp(-t*drag_coeff)
def main():
"""
Run test where the trajectory of initially moving particles in a
quiescent flow is compared to the analytical exponential
Run test where the trajectory of initially moving particles with varying drag coefficient in a
quiescent flow is compared to the analytically expected exponential. The generated plot shows the accuracy of the corresponding particles for a certain fixed timestep.
"""
drag_coefficients = ['30.','20.','10.','1.','0.1']
drag_coefficients = ['30','20','10','1']
niterations = 100
nparticles = 1
njobs = 1
dt=0.01
fig = plt.figure()
plt.ylabel('$v$')
plt.xlabel('$t$')
plt.yscale('log')
j=0
for drag_coeff in drag_coefficients:
c = DNS()
......@@ -68,7 +66,6 @@ def main():
'--niter_todo', '{0}'.format(niterations),
'--niter_out', '{0}'.format(1),
'--niter_stat', '1',
'--checkpoints_per_file', '{0}'.format(3),
'--nparticles', '{0}'.format(nparticles),
'--initial_field_amplitude', '0',
'--initial_particle_vel', '0.05',
......@@ -77,29 +74,25 @@ def main():
'--njobs', '{0}'.format(njobs),
'--wd', './'] +
sys.argv[1:])
f = h5py.File('quiescent_nsve_stokes_particles_drag'+str(j)+'_particles.h5', 'r')
f2 = h5py.File('quiescent_nsve_stokes_particles_drag'+str(j)+'.h5', 'r')
particle_momentum = [np.sqrt(f['tracers0/momentum/'+'{}'.format(i)][0][0]**2
+f['tracers0/momentum/'+'{}'.format(i)][0][1]**2
+f['tracers0/momentum/'+'{}'.format(i)][0][2]**2)
for i in range(niterations)]
time_values = np.arange(0, niterations*dt, dt)
def theory_exponential(t):
return f2['parameters/initial_particle_vel']*np.exp(-t*np.float(drag_coeff))
numerics, = plt.plot(time_values, particle_momentum, alpha = 0.3, lw=2)
plt.plot(time_values, theory_exponential(time_values), color=numerics.get_color(), ls = '--')
plt.plot(0,0,color=numerics.get_color(),label=r'$\frac{1}{St}=$'+drag_coeff)
plt.legend(loc='lower left')
print(f['tracers0/momentum/0'][:])
print(f['tracers0/momentum/100'][:]-f['tracers0/momentum/0'][:])
print(f['tracers0/position/100'][:]-f['tracers0/position/0'][:])
plt.savefig('Stokes_test1.pdf')
numerics, = plt.plot(time_values, particle_momentum, alpha = 1, lw=2)
plt.plot(time_values,
theory_exponential(time_values, f2['parameters/initial_particle_vel'],np.float(drag_coeff)),
color='black', ls = '--')
plt.plot(0,0,color=numerics.get_color(),label=r'drag coefficient '+drag_coeff)
j+=1
plt.plot(0,
0,
color='black', ls = '--',
label='analytic sol.')
plt.legend(loc='lower left')
plt.savefig('Stokes_quiescent_test.pdf')
return None
if __name__ == '__main__':
......
......@@ -46,8 +46,6 @@ int NSVE_Stokes_particles<rnumber>::initialize(void)
particles_inner_computer_2nd_order_Stokes<double, long long int> current_particles_inner_computer;
current_particles_inner_computer.set_drag_coefficient(this->drag_coefficient);
//DEBUG_MSG("drag coefficient is set to %f \n", current_particles_inner_computer.get_drag_coefficient());
//DEBUG_MSG_WAIT(MPI_COMM_WORLD, "before call to particles_system_builder\n");
this->ps = particles_system_builder_with_p2p(
this->fs->cvelocity, // (field object)
this->fs->kk, // (kspace object, contains dkx, dky, dkz)
......@@ -79,7 +77,6 @@ int NSVE_Stokes_particles<rnumber>::initialize(void)
"tracers0",
"position/0");
this->particles_sample_writer_mpi->setParticleFileLayout(this->ps->getParticleFileLayout());
//DEBUG_MSG("drag coefficient is after initialization %f \n", current_particles_inner_computer.get_drag_coefficient());
return EXIT_SUCCESS;
}
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment