 |
tdouble_gyre.jl - Granular.jl - Julia package for granular dynamics simulation |
 |
git clone git://src.adamsgaard.dk/Granular.jl |
 |
Log |
|
Files |
|
Refs |
|
README |
|
LICENSE |
|
--- |
|
tdouble_gyre.jl (3099B) |
|
--- |
|
1 #!/usr/bin/env julia |
|
2 import Granular |
|
3 using Random |
|
4 |
|
5 sim = Granular.createSimulation(id="double_gyre") |
|
6 |
|
7 # Initialize ocean |
|
8 L = [100e3, 50e3, 1e3] |
|
9 Ly_constriction = 20e3 |
|
10 #n = [750, 500, 2] # high resolution |
|
11 n = [30, 15, 2] # intermedite resolution |
|
12 #n = [8, 5, 2] # coarse resolution |
|
13 sim.ocean = Granular.createRegularOceanGrid(n, L, name="double_gyre") |
|
14 |
|
15 epsilon = 0.25 # amplitude of periodic oscillations |
|
16 t = 0. |
|
17 a = epsilon*sin(2.0*pi*t) |
|
18 b = 1.0 - 2.0*epsilon*sin(2.0*pi*t) |
|
19 for i=1:size(sim.ocean.u, 1) |
|
20 for j=1:size(sim.ocean.u, 2) |
|
21 |
|
22 x = sim.ocean.xq[i, j]/(L[1]*0.5) # x in [0;2] |
|
23 y = sim.ocean.yq[i, j]/L[2] # y in [0;1] |
|
24 |
|
25 f = a*x^2.0 + b*x |
|
26 df_dx = 2.0*a*x + b |
|
27 |
|
28 sim.ocean.u[i, j, 1, 1] = -pi/10.0*sin(pi*f)*cos(pi*y) * 1e1 |
|
29 sim.ocean.v[i, j, 1, 1] = pi/10.0*cos(pi*f)*sin(pi*y)*df_dx * 1e1 |
|
30 end |
|
31 end |
|
32 |
|
33 # Initialize confining walls, which are ice floes that are fixed in space |
|
34 r = minimum(L[1:2]./n[1:2])/2.0 |
|
35 h = 1. |
|
36 |
|
37 ## N-S wall segments |
|
38 for y in range(r, stop=L[2]-r, length=Int(round((L[2] - 2.0*r)/(r*2)))) |
|
39 Granular.addGrainCylindrical!(sim, [r, y], r, h, fixed=true, |
|
40 verbose=false) |
|
41 Granular.addGrainCylindrical!(sim, [L[1]-r, y], r, h, fixed=true, |
|
42 verbose=false) |
|
43 end |
|
44 |
|
45 ## E-W wall segments |
|
46 for x in range(3.0*r, stop=L[1]-3.0*r, |
|
47 length=Int(round((L[1] - 6.0*r)/(r*2)))) |
|
48 Granular.addGrainCylindrical!(sim, [x, r], r, h, fixed=true, |
|
49 verbose=false) |
|
50 Granular.addGrainCylindrical!(sim, [x, L[2]-r], r, h, fixed=true, |
|
51 verbose=false) |
|
52 end |
|
53 |
|
54 n_walls = length(sim.grains) |
|
55 @info "added $(n_walls) fixed ice floes as walls" |
|
56 |
|
57 |
|
58 |
|
59 # Initialize ice floes everywhere |
|
60 floe_padding = 0.5*r |
|
61 noise_amplitude = 0.8*floe_padding |
|
62 Random.seed!(1) |
|
63 for y in (4.0*r + noise_amplitude):(2.0*r + floe_padding):(L[2] - 4.0*r … |
|
64 noise_amplitu… |
|
65 |
|
66 for x in (4.0*r + noise_amplitude):(2.0*r + floe_padding):(L[1] - 4.… |
|
67 noise_amp… |
|
68 #if iy % 2 == 0 |
|
69 #x += 1.5*r |
|
70 #end |
|
71 |
|
72 x_ = x + noise_amplitude*(0.5 - rand()) |
|
73 y_ = y + noise_amplitude*(0.5 - rand()) |
|
74 |
|
75 Granular.addGrainCylindrical!(sim, [x_, y_], r, h, verbose=false) |
|
76 end |
|
77 end |
|
78 n = length(sim.grains) - n_walls |
|
79 @info "added $n ice floes" |
|
80 |
|
81 # Remove old simulation files |
|
82 Granular.removeSimulationFiles(sim) |
|
83 |
|
84 k_n = 1e6 # N/m |
|
85 gamma_t = 1e7 # N/(m/s) |
|
86 mu_d = 0.7 |
|
87 rotating = false |
|
88 for i=1:length(sim.grains) |
|
89 sim.grains[i].contact_stiffness_normal = k_n |
|
90 sim.grains[i].contact_stiffness_tangential = k_n |
|
91 sim.grains[i].contact_viscosity_tangential = gamma_t |
|
92 sim.grains[i].contact_dynamic_friction = mu_d |
|
93 sim.grains[i].rotating = rotating |
|
94 end |
|
95 |
|
96 # Set temporal parameters |
|
97 Granular.setTotalTime!(sim, 12.0*60.0*60.0) |
|
98 Granular.setOutputFileInterval!(sim, 60.0) |
|
99 Granular.setTimeStep!(sim) |
|
100 |
|
101 Granular.run!(sim) |
