atmosphere.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)