|
59 | 59 | q_prim_vf(advxb)%sf(i, j, k) = patch_icpp(1)%alpha(1) |
60 | 60 | q_prim_vf(advxe)%sf(i, j, k) = patch_icpp(1)%alpha(2) |
61 | 61 | end if |
62 | | - |
| 62 | + |
63 | 63 | case (302) ! (3D lung geometry in X direction - axisym, with smoothing) |
64 | | - lam = 200.E-06_wp |
| 64 | + lam = 200.e-06_wp |
65 | 65 | amp = patch_icpp(patch_id)%a2 |
66 | 66 | h = 0.125_wp*amp |
67 | | - |
| 67 | + |
68 | 68 | intH = amp/2._wp*(sin(2._wp*pi*y_cc(j)/lam + pi/2._wp) + sin(2._wp*pi*z_cc(k)/lam + pi/2._wp)) |
69 | | - |
70 | | - alph = patch_icpp(2)%alpha(1) + (patch_icpp(1)%alpha(1)-patch_icpp(2)%alpha(1))/(h)*(x_cc(i)-(intH-h/2._wp)) |
71 | | - |
| 69 | + |
| 70 | + alph = patch_icpp(2)%alpha(1) + (patch_icpp(1)%alpha(1) - patch_icpp(2)%alpha(1))/(h)*(x_cc(i) - (intH - h/2._wp)) |
| 71 | + |
72 | 72 | if (x_cc(i) > intH + h/2) then |
73 | | - |
| 73 | + |
74 | 74 | q_prim_vf(advxb)%sf(i, j, k) = patch_icpp(1)%alpha(1) |
75 | 75 | q_prim_vf(advxe)%sf(i, j, k) = patch_icpp(1)%alpha(2) |
76 | 76 | q_prim_vf(contxb)%sf(i, j, k) = patch_icpp(1)%alpha_rho(1) |
77 | 77 | q_prim_vf(contxe)%sf(i, j, k) = patch_icpp(1)%alpha_rho(2) |
78 | 78 | q_prim_vf(E_idx)%sf(i, j, k) = patch_icpp(1)%pres |
79 | | - |
80 | | - else if ((x_cc(i) .le. intH + h/2) .and. (x_cc(i) .ge. intH - h/2._wp)) then |
81 | | - |
| 79 | + |
| 80 | + else if ((x_cc(i) <= intH + h/2) .and. (x_cc(i) >= intH - h/2._wp)) then |
| 81 | + |
82 | 82 | q_prim_vf(advxb)%sf(i, j, k) = alph !0.5 |
83 | 83 | q_prim_vf(advxe)%sf(i, j, k) = 1._wp - alph !0.5 |
84 | 84 | q_prim_vf(contxb)%sf(i, j, k) = patch_icpp(1)%alpha_rho(1)/patch_icpp(1)%alpha(1)*alph!0.5 |
85 | | - q_prim_vf(contxe)%sf(i, j, k) = patch_icpp(2)%alpha_rho(2)/patch_icpp(2)%alpha(2)*(1-alph)!0.5 |
| 85 | + q_prim_vf(contxe)%sf(i, j, k) = patch_icpp(2)%alpha_rho(2)/patch_icpp(2)%alpha(2)*(1 - alph)!0.5 |
86 | 86 | q_prim_vf(E_idx)%sf(i, j, k) = patch_icpp(1)%pres |
87 | | - |
| 87 | + |
88 | 88 | end if |
89 | | - |
| 89 | + |
90 | 90 | case (303) ! pre_stress for hyperelasticity, bubble in material |
91 | 91 |
|
92 | | - R0ref = 30E-6_wp ! equilibrium radius |
| 92 | + R0ref = 30e-6_wp ! equilibrium radius |
93 | 93 | Rinit = patch_icpp(3)%radius ! initial radius |
94 | 94 | x_bcen = patch_icpp(3)%x_centroid |
95 | 95 | y_bcen = patch_icpp(3)%y_centroid |
|
101 | 101 | phi = atan2(y_ccs, x_ccs) |
102 | 102 | theta = atan2(sqrt(x_ccs**2._wp + y_ccs**2._wp), z_ccs) |
103 | 103 | !spherical coord, assuming Rmax=1 |
104 | | - xi_sph = (rcoord**3._wp - R0ref**3._wp + Rinit**3._wp)**(1._wp/3._wp) |
| 104 | + xi_sph = (rcoord**3._wp - R0ref**3._wp + Rinit**3._wp)**(1._wp/3._wp) |
105 | 105 | xi_cart(1) = xi_sph*sin(theta)*cos(phi) |
106 | 106 | xi_cart(2) = xi_sph*sin(theta)*sin(phi) |
107 | 107 | xi_cart(3) = xi_sph*cos(theta) |
|
0 commit comments