adds initialisation for Stokes particles

TurTLE/DNS.py

@@ -703,6 +703,9 @@ class DNS(_code):
             pars['spectrum_slope'] = float(-5./3)
             pars['spectrum_k_cutoff'] = float(16)
             pars['spectrum_coefficient'] = float(1)
+        if dns_type == 'NSVE_Stokes_particles':
+            pars['initial_field_amplitude'] = float(0.05)
+            pars['initial_particle_vel'] = float(0.05)
         return pars
     def prepare_launch(
             self,
@@ -897,11 +900,17 @@ class DNS(_code):
                     if nn > batch_size:
                         dset[cc*batch_size:(cc+1)*batch_size, :3] = get_random_phases(batch_size)
                         if dset.shape[1] == 6:
+                            if self.dns_type == 'NSVE_Stokes_particles':
+                                dset[cc*batch_size:(cc+1)*batch_size, 3:] = self.parameters['initial_particle_vel']*get_random_versors(batch_size)
+                            else:
                                 dset[cc*batch_size:(cc+1)*batch_size, 3:] = get_random_versors(batch_size)                             
                         nn -= batch_size
                     else:
                         dset[cc*batch_size:cc*batch_size+nn, :3] = get_random_phases(nn)
                         if dset.shape[1] == 6:
+                            if self.dns_type == 'NSVE_Stokes_particles':
+                                dset[cc*batch_size:cc*batch_size+nn, 3:] =  self.parameters['initial_particle_vel']*get_random_versors(batch_size)
+                            else:
                                 dset[cc*batch_size:cc*batch_size+nn, 3:] = get_random_versors(nn)
                         nn = 0
                     cc += 1
@@ -1106,6 +1115,12 @@ class DNS(_code):
                         'vorticity/complex/{0}'.format(opt.src_iteration),
                         f,
                         'vorticity/complex/{0}'.format(0))
+            else:
+                if self.dns_type == 'NSVE_Stokes_particles':
+                  data = self.generate_vector_field(
+                       write_to_file = False,
+                       spectra_slope = 2.0,
+                       amplitude = self.parameters['initial_field_amplitude'])
                 else:
                     data = self.generate_vector_field(
                        write_to_file = False,