diff --git a/examples/sediment_deposit/particles_above_salt_bc.py b/examples/sediment_deposit/sediment_deposit.py similarity index 94% rename from examples/sediment_deposit/particles_above_salt_bc.py rename to examples/sediment_deposit/sediment_deposit.py index d56d9d786bc641df8bcd6c406a809e24b2332b41..54b493b9d1aefec28ec4580ef68199697f998a9f 100644 --- a/examples/sediment_deposit/particles_above_salt_bc.py +++ b/examples/sediment_deposit/sediment_deposit.py @@ -17,8 +17,17 @@ def init_velocity(data, coords, component=None): for d in data: d[...] = 0.0 -def init_sediment(data, coords, l0): - raise NotImplementedError +def init_sediment(data, coords, nblobs, rblob): + r2 = rblob * rblob + X = 7*np.random.rand(nblobs) + Y = 1*np.random.rand(nblobs) + data = data[0] + coords = coords[0] + + data[...] = 0.0 + for (x,y) in zip(X,Y): + D = (coords[0]-x)**2 + (coords[1]-y)**2 + data[D<r2] = 0.5 def compute(args): from hysop import Field, Box, Simulation, Problem, MPIParams, IOParams, vprint, \ @@ -53,22 +62,18 @@ def compute(args): # Constants dim = args.ndim if (dim==2): - Sc = 0.70 Xo = (0.0,)*dim - Xn = (1.0,)*dim - N = args.npts + Xn = (1.0,7.0) + nblobs = 240 + rblob = 1e-2 + npts = args.npts else: msg='The {}D has not been implemented yet.'.format(dim) raise NotImplementedError(msg) - nu_S = ScalarParameter(name='nu_S', dtype=args.dtype, const=True, initial_value=0.01) - nu_W = ScalarParameter(name='nu_W', dtype=args.dtype, const=True, initial_value=0.01) + nu_S = ScalarParameter(name='nu_S', dtype=args.dtype, const=True, initial_value=1.00) + nu_W = ScalarParameter(name='nu_W', dtype=args.dtype, const=True, initial_value=1.00) - # Define the domain - #npts=N[::-1] - #Xo=Xo[::-1] - #Xn=Xn[::-1] - lboundaries = (BoxBoundaryCondition.PERIODIC,)*(dim-1)+(BoxBoundaryCondition.SYMMETRIC,) rboundaries = (BoxBoundaryCondition.PERIODIC,)*(dim-1)+(BoxBoundaryCondition.SYMMETRIC,) @@ -218,7 +223,7 @@ def compute(args): ### Adaptive timestep operator adapt_dt = AdaptiveTimeStep(dt, equivalent_CFL=True, - name='merge_dt', pretty_name='dt', ) + name='merge_dt', pretty_name='dt', max_dt=0.1) dt_cfl = adapt_dt.push_cfl_criteria(cfl=args.cfl, Fmin=min_max_U.Fmin, Fmax=min_max_U.Fmax, @@ -272,7 +277,7 @@ def compute(args): # Initialize vorticity, velocity, S on all topologies problem.initialize_field(field=velo, formula=init_velocity) problem.initialize_field(field=vorti, formula=init_vorticity) - problem.initialize_field(field=S, formula=init_sediment) + problem.initialize_field(field=S, formula=init_sediment, nblobs=nblobs, rblob=rblob) # Finally solve the problem problem.solve(simu, dry_run=args.dry_run) @@ -310,12 +315,12 @@ if __name__=='__main__': parser = ParticleAboveSaltArgParser() - parser.set_defaults(impl='cl', ndim=2, npts=(500,3000), + parser.set_defaults(impl='cl', ndim=2, npts=(500,3001), box_origin=(0.0,), box_length=(1.0,), tstart=0.0, tend=100.0, dt=1e-6, cfl=4.00, lcfl=0.95, - dump_times=tuple(float(x) for x in range(0,500,10)), - dump_freq=0) + dump_times=tuple(float(x) for x in range(0,100,10)), + dump_freq=1) parser.run(compute)