tatmosphere.jl - Granular.jl - Julia package for granular dynamics simulation
(HTM) git clone git://src.adamsgaard.dk/Granular.jl
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) README
(DIR) LICENSE
---
tatmosphere.jl (6765B)
---
1 #!/usr/bin/env julia
2
3 # Check if atmosphere-specific functions and grid operations are functioning
4 # correctly
5
6 @info "Testing regular grid generation"
7 sim = Granular.createSimulation()
8 sim.atmosphere = Granular.createRegularAtmosphereGrid([6, 6, 6], [1., 1., 1.])
9 @test size(sim.atmosphere.xq) == (7, 7)
10 @test size(sim.atmosphere.yq) == (7, 7)
11 @test size(sim.atmosphere.xh) == (6, 6)
12 @test size(sim.atmosphere.yh) == (6, 6)
13 @test sim.atmosphere.xq[1, :, 1] ≈ zeros(7)
14 @test sim.atmosphere.xq[4, :, 1] ≈ .5 * ones(7)
15 @test sim.atmosphere.xq[end, :, 1] ≈ 1. * ones(7)
16 @test sim.atmosphere.yq[:, 1, 1] ≈ zeros(7)
17 @test sim.atmosphere.yq[:, 4, 1] ≈ .5 * ones(7)
18 @test sim.atmosphere.yq[:, end, 1] ≈ 1. * ones(7)
19 @test size(sim.atmosphere.u) == (7, 7, 6, 1)
20 @test size(sim.atmosphere.v) == (7, 7, 6, 1)
21 @test sim.atmosphere.u ≈ zeros(7, 7, 6, 1)
22 @test sim.atmosphere.v ≈ zeros(7, 7, 6, 1)
23
24 @info "Testing velocity drag interaction (static atmosphere)"
25 Granular.addGrainCylindrical!(sim, [.5, .5], .25, .1)
26 Granular.setTotalTime!(sim, 5.)
27 Granular.setTimeStep!(sim)
28 sim_init = deepcopy(sim)
29 sim.grains[1].lin_vel[1] = 0.1
30 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
31 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
32 Granular.run!(sim, verbose=false)
33 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
34 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
35 @test E_kin_rot_init ≈ E_kin_rot_final # no rotation before or after
36 @test E_kin_lin_init > E_kin_lin_final # linear velocity lost due to atmos drag
37 @test sim.grains[1].atmosphere_stress[1] < 0.
38 @test sim.grains[1].atmosphere_stress[2] ≈ 0.
39
40 @info "Testing velocity drag interaction (static ice floe)"
41 sim = deepcopy(sim_init)
42 sim.atmosphere.v[:, :, 1, 1] .= 0.1
43 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
44 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
45 Granular.run!(sim, verbose=false)
46 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
47 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
48 @test E_kin_rot_init ≈ E_kin_rot_final # no rotation before or after
49 @test E_kin_lin_init < E_kin_lin_final # linear vel. gained due to atmos drag
50 @test sim.grains[1].atmosphere_stress[1] ≈ 0.
51 @test sim.grains[1].atmosphere_stress[2] > 0.
52
53 @info "Testing vortex interaction (static atmosphere)"
54 sim = deepcopy(sim_init)
55 sim.grains[1].ang_vel[3] = 0.1
56 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
57 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
58 Granular.run!(sim, verbose=false)
59 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
60 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
61 @test E_kin_rot_init > E_kin_rot_final # energy lost to atmosphere
62 @test sim.grains[1].ang_vel[3] > 0. # check angular velocity orientation
63 @test sim.grains[1].ang_pos[3] > 0. # check angular position orientation
64 @test E_kin_lin_init ≈ E_kin_lin_final # no linear velocity gained
65
66 @info "Testing vortex interaction (static ice floe)"
67 sim = deepcopy(sim_init)
68 sim.atmosphere = Granular.createRegularAtmosphereGrid([1, 1, 1], [1., 1., 1.])
69 sim.grains[1].lin_pos[1] = 0.5
70 sim.grains[1].lin_pos[2] = 0.5
71 sim.atmosphere.v[1, :, 1, 1] .= -0.1
72 sim.atmosphere.v[2, :, 1, 1] .= 0.1
73 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
74 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
75 Granular.run!(sim, verbose=false)
76 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
77 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
78 @test sim.grains[1].ang_vel[3] > 0. # check angular velocity orientation
79 @test sim.grains[1].ang_pos[3] > 0. # check angular position orientation
80 @test E_kin_rot_init < E_kin_rot_final # rotation after due to atm vortex
81 @test E_kin_lin_init ≈ E_kin_lin_final # no linear velocity gained
82
83 sim = deepcopy(sim_init)
84 sim.atmosphere = Granular.createRegularAtmosphereGrid([1, 1, 1], [1., 1., 1.])
85 sim.grains[1].lin_pos[1] = 0.5
86 sim.grains[1].lin_pos[2] = 0.5
87 sim.atmosphere.v[1, :, 1, 1] .= 0.1
88 sim.atmosphere.v[2, :, 1, 1] .= -0.1
89 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
90 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
91 Granular.run!(sim, verbose=false)
92 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
93 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
94 @test sim.grains[1].ang_vel[3] < 0. # check angular velocity orientation
95 @test sim.grains[1].ang_pos[3] < 0. # check angular position orientation
96 @test E_kin_rot_init < E_kin_rot_final # rotation after due to atm vortex
97 @test E_kin_lin_init ≈ E_kin_lin_final # no linear velocity gained
98
99 sim = deepcopy(sim_init)
100 sim.atmosphere = Granular.createRegularAtmosphereGrid([1, 1, 1], [1., 1., 1.])
101 sim.grains[1].lin_pos[1] = 0.5
102 sim.grains[1].lin_pos[2] = 0.5
103 sim.atmosphere.u[:, 1, 1, 1] .= -0.1
104 sim.atmosphere.u[:, 2, 1, 1] .= 0.1
105 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
106 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
107 Granular.run!(sim, verbose=false)
108 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
109 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
110 @test sim.grains[1].ang_vel[3] < 0. # check angular velocity orientation
111 @test sim.grains[1].ang_pos[3] < 0. # check angular position orientation
112 @test E_kin_rot_init < E_kin_rot_final # rotation after due to atm vortex
113 @test E_kin_lin_init ≈ E_kin_lin_final # no linear velocity gained
114
115 sim = deepcopy(sim_init)
116 sim.atmosphere = Granular.createRegularAtmosphereGrid([1, 1, 1], [1., 1., 1.])
117 sim.grains[1].lin_pos[1] = 0.5
118 sim.grains[1].lin_pos[2] = 0.5
119 sim.atmosphere.u[:, 1, 1, 1] .= 0.1
120 sim.atmosphere.u[:, 2, 1, 1] .= -0.1
121 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
122 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
123 Granular.run!(sim, verbose=false)
124 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
125 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
126 @test sim.grains[1].ang_vel[3] > 0. # check angular velocity orientation
127 @test sim.grains[1].ang_pos[3] > 0. # check angular position orientation
128 @test E_kin_rot_init < E_kin_rot_final # rotation after due to atm vortex
129 @test E_kin_lin_init ≈ E_kin_lin_final # no linear velocity gained
130
131 sim = Granular.createSimulation()
132 sim.atmosphere = Granular.createRegularAtmosphereGrid([6, 6, 6], [1., 1., 1.])
133 sim2 = Granular.createSimulation()
134 sim2.atmosphere = Granular.createRegularAtmosphereGrid([6, 6, 6], [1., 1., 1.])
135 Granular.compareAtmospheres(sim.atmosphere, sim2.atmosphere)