41 un_m1, un, pn_m1, pn, phin_m1, phin)
44 REAL(KIND=8),
INTENT(OUT):: time
45 REAL(KIND=8),
INTENT(IN) :: dt
46 INTEGER,
DIMENSION(:),
INTENT(IN) :: list_mode
47 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(OUT):: un_m1, un
48 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(OUT):: pn_m1, pn, phin_m1, phin
50 REAL(KIND=8),
DIMENSION(mesh_c%np) :: pn_m2
53 DO i= 1,
SIZE(list_mode)
57 un_m1(:,j,i) =
vv_exact(j,mesh_f%rr,mode,time-dt)
58 un(:,j,i) =
vv_exact(j,mesh_f%rr,mode,time)
62 pn_m2(:) =
pp_exact(j,mesh_c%rr,mode,time-2*dt)
63 pn_m1(:,j,i) =
pp_exact(j,mesh_c%rr,mode,time-dt)
64 pn(:,j,i) =
pp_exact(j,mesh_c%rr,mode,time)
65 phin_m1(:,j,i) = pn_m1(:,j,i) - pn_m2(:)
66 phin(:,j,i) = pn(:,j,i) - pn_m1(:,j,i)
75 REAL(KIND=8),
INTENT(OUT):: time
76 REAL(KIND=8),
INTENT(IN) :: dt
77 INTEGER,
DIMENSION(:),
INTENT(IN) :: list_mode
78 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(OUT):: tempn_m1, tempn
82 DO i= 1,
SIZE(list_mode)
117 INTEGER ,
INTENT(IN) :: TYPE
118 REAL(KIND=8),
DIMENSION(:,:),
INTENT(IN) :: rr
119 INTEGER ,
INTENT(IN) :: mode, i
120 REAL(KIND=8),
INTENT(IN) :: time
121 REAL(KIND=8),
INTENT(IN) :: Re
122 CHARACTER(LEN=2),
INTENT(IN) :: ty
123 REAL(KIND=8),
DIMENSION(:,:,:),
OPTIONAL,
INTENT(IN) :: opt_density
124 REAL(KIND=8),
DIMENSION(:,:,:),
OPTIONAL,
INTENT(IN) :: opt_tempn
125 REAL(KIND=8),
DIMENSION(SIZE(rr,2)) :: vv, r, z
126 REAL(KIND=8) :: alpha, r0 = 0.5d0
127 CHARACTER(LEN=2) :: np
129 IF (
PRESENT(opt_density))
CALL error_petsc(
'density should not be present for test 30')
131 alpha =
inputs%gravity_coefficient
137 ELSE IF (type==6)
THEN 145 vv = vv -(2*r*(r**4 - 2*r**3*r0 + r0**2 + r**2*(-3 + r0**2))*cos(time)*cos(z) + &
146 (r - r0)*re*cos(time)**2*(14*r**3 - 5*r**2*r0 - 5*r*r0**2 + 2*r0**3 + &
147 (9*r**5 - 21*r**4*r0 + 5*r*r0**2 - 2*r0**3 + r**3*(-14 + 15*r0**2) + r**2*(5*r0 - 3*r0**3))*cos(2*z)) - &
148 2*r**3*(r - r0)**2*re*cos(z)*sin(time))/(2.*r**3*re)
149 ELSE IF (mode==1)
THEN 150 vv = vv +(-(r*(r**4 - 2*r*r0 - 2*r**3*r0 + 2*r0**2 + r**2*(-2 + r0**2))*cos(time)*cos(z)) - &
151 (3*r**2 - 4*r*r0 + r0**2)*re*cos(time)**2*(3*r**2 - r0**2 + (2*r**4 - 4*r**3*r0 + r0**2 + &
152 r**2*(-3 + 2*r0**2))*cos(2*z)) + r**3*(r - r0)**2*re*cos(z)*sin(time))/ (r**3*re)
153 ELSE IF (mode==2)
THEN 154 vv = vv -((r - r0)*cos(time)**2*(4*r**2 - r*r0 - r0**2 + (3*r**4 - 7*r**3*r0 + r0**2 + &
155 r**2*(-4 + 5*r0**2) + r*(r0 - r0**3))*cos(2*z)))/(2.*r**2)
157 ELSE IF (type==2)
THEN 159 vv = vv + ((r - r0)**2*cos(time)*(-2*r*cos(z) + re*cos(time)*(r**4 - r*r0 - 2*r**3*r0 + &
160 r0**2 + r**2*(-3 + r0**2) + (3*r**2 + r*r0 - r0**2)*cos(2*z))))/(r**3*re)
161 ELSE IF (mode==2)
THEN 162 vv = vv + ((r - r0)**2*cos(time)**2*(r**4 - r*r0 - 2*r**3*r0 + r0**2 + r**2*(-3 + r0**2) + &
163 (3*r**2 + r*r0 - r0**2)*cos(2*z)))/(2.*r**3)
165 ELSE IF (type==3)
THEN 167 vv = vv + (-3*r*(r - r0)**3*(3*r - r0)*re*cos(time)**2 + 2*(r**4 - 2*r**3*r0 + r0**2 + &
168 r**2*(-3 + r0**2))*cos(time)*cos(z) - 2*r**2*(r - r0)**2*re*cos(z)*sin(time))/(2.*r**2*re)
169 ELSE IF (mode==1)
THEN 170 vv = vv -(-2*r*(r**4 - 2*r*r0 - 2*r**3*r0 + 2*r0**2 + r**2*(-2 + r0**2))*cos(time)*cos(z) + &
171 (r - r0)**3*(3*r - r0)*re*cos(time)**2*(1 + 4*r**2 - cos(2*z)) + &
172 2*r**3*(r - r0)**2*re*cos(z)*sin(time))/(2.*r**3*re)
173 ELSE IF (mode==2)
THEN 174 vv = vv + ((r - r0)**3*(3*r - r0)*cos(time)**2*(-1 - r**2 + cos(2*z)))/(2.*r**3)
176 ELSE IF (type==4)
THEN 178 vv = vv + ((r - r0)**2*cos(time)*(-4*r*cos(z) + re*cos(time)*((-3*r + r0)**2 + &
179 (2*r**4 + 6*r*r0 - 4*r**3*r0 - r0**2 + r**2*(-9 + 2*r0**2))*cos(2*z))))/(2.*r**3*re)
180 ELSE IF (mode==2)
THEN 181 vv = vv + ((r - r0)**2*cos(time)**2*(4*r - 2*r0 + (r*(-4 + (r - r0)**2) + 2*r0)*cos(2*z)))/(2.*r**2)
183 ELSE IF (type==5)
THEN 185 vv = vv + (((3*r**4 - 4*r**3*r0 - r0**2 + r**2*(-3 + r0**2))*cos(time) + &
186 r*(r - r0)**2*(3*r**4 - r*r0 - 6*r**3*r0 + 2*r0**2 + r**2*(-4 + 3*r0**2))*re*cos(time)**2*cos(z) - &
187 r**2*(3*r**2 - 4*r*r0 + r0**2)*re*sin(time))*sin(z))/(r**3*re)
188 ELSE IF (mode==1)
THEN 189 vv = vv + ((3*r**3 - 4*r0 - 4*r**2*r0 + r*r0**2)*cos(time)*sin(z) - r*(3*r**2 - 4*r*r0 + r0**2)*re*sin(time)*sin(z) + &
190 ((r - r0)**2*(4*r**4 - 2*r*r0 - 8*r**3*r0 + 3*r0**2 + r**2*(-5 + 4*r0**2))*re*cos(time)**2*sin(2*z))/2.)/(r**2*re)
191 ELSE IF (mode==2)
THEN 192 vv = vv + ((r - r0)**2*(r**4 - r*r0 - 2*r**3*r0 + r0**2 + r**2*(-1 + r0**2))*cos(time)**2*sin(2*z))/(2.*r**2)
194 ELSE IF (type==6)
THEN 196 vv = vv -((((-r**3 + 2*r0 + 2*r**2*r0 - r*r0**2)*cos(time) + 4*r*(r - r0)**3*re*cos(time)**2*cos(z) + &
197 r*(r - r0)**2*re*sin(time))*sin(z))/(r**2*re))
198 ELSE IF (mode==2)
THEN 199 vv = vv -(((r - r0)**3*cos(time)**2*sin(2*z))/r)
213 INTEGER ,
INTENT(IN) :: TYPE
214 REAL(KIND=8),
DIMENSION(:,:),
INTENT(IN) :: rr
215 INTEGER ,
INTENT(IN) :: m
216 REAL(KIND=8),
INTENT(IN) :: t
217 REAL(KIND=8),
DIMENSION(SIZE(rr,2)) :: vv, r, z, c, lambda
219 REAL(KIND=8) :: r0 = 0.5d0
226 c(i) =
inputs%vol_heat_capacity(1)
228 c(i) =
inputs%vol_heat_capacity(2)
234 lambda(i) =
inputs%temperature_diffusivity(1)
236 lambda(i) =
inputs%temperature_diffusivity(2)
242 vv = ((-9*r + r**3 + 4*r0 - r**2*r0) * lambda * cos(t) + &
243 c * r**2 * (-r + r0) * sin(t)) * sin(z) / lambda
245 vv = ((-8*r + r**3 + 3*r0 - r**2*r0) * lambda * cos(t) + &
246 c * r**2 * (-r + r0) * sin(t)) * sin(z) / lambda
258 vv(i) = vv(i) + (3*c(i) * r(i) * (r(i) - r0)**2 * r0 &
259 * cos(t)**2 * sin(2*z(i))) / (4 * lambda(i))
265 vv(i) = vv(i) + (c(i) * r(i) * (r(i) - r0)**2 * r0 &
266 * cos(t)**2 * sin(2*z(i))) / lambda(i)
272 vv(i) = vv(i) + (c(i) * r(i) * (r(i) - r0)**2 * r0 &
273 * cos(t)**2 * sin(2*z(i))) / (4 *lambda(i))
297 FUNCTION vv_exact(TYPE,rr,m,t) RESULT(vv)
299 INTEGER ,
INTENT(IN) :: TYPE
300 REAL(KIND=8),
DIMENSION(:,:),
INTENT(IN) :: rr
301 INTEGER,
INTENT(IN) :: m
302 REAL(KIND=8),
INTENT(IN) :: t
303 REAL(KIND=8),
DIMENSION(SIZE(rr,2)) :: vv, r, z
304 REAL(KIND=8) :: r0 = 0.5d0
310 IF ((m==0).OR.(m==1))
THEN 311 vv = -(r-r0)**2*cos(t)*cos(z)
315 ELSE IF (type==3)
THEN 316 IF ((m==0).OR.(m==1))
THEN 317 vv = (r-r0)**2*cos(t)*cos(z)
321 ELSE IF (type==5)
THEN 322 IF ((m==0).OR.(m==1))
THEN 323 vv = (r-r0)*cos(t)*sin(z)/r * (3*r-r0)
327 ELSE IF (type==6)
THEN 329 vv = (r-r0)*cos(t)*sin(z)/r * (r-r0)
357 FUNCTION pp_exact(TYPE,rr,m,t) RESULT (vv)
359 INTEGER ,
INTENT(IN) :: TYPE
360 REAL(KIND=8),
DIMENSION(:,:),
INTENT(IN) :: rr
361 INTEGER ,
INTENT(IN) :: m
362 REAL(KIND=8),
INTENT(IN) :: t
363 REAL(KIND=8),
DIMENSION(SIZE(rr,2)) :: vv
371 n=type; n=
SIZE(rr,1); n=m; r=t
378 INTEGER ,
INTENT(IN) :: TYPE
379 REAL(KIND=8),
DIMENSION(:,:),
INTENT(IN) :: rr
380 INTEGER ,
INTENT(IN) :: m
381 REAL(KIND=8),
INTENT(IN) :: t
382 REAL(KIND=8),
DIMENSION(SIZE(rr,2)) :: vv, r, z, lambda
383 REAL(KIND=8) :: r0=0.5d0
391 lambda(i) =
inputs%temperature_diffusivity(1)
393 lambda(i) =
inputs%temperature_diffusivity(2)
398 IF ((type==1).AND.((m==0).OR.(m==1)))
THEN 399 vv = r**2*(r-r0)*sin(z)*cos(t) / lambda
446 INTEGER ,
INTENT(IN) ::
TYPE, n_start
447 INTEGER,
INTENT(IN) :: mode
448 REAL(KIND=8),
INTENT(IN) :: t
449 REAL(KIND=8),
DIMENSION(H_Mesh%np) :: vv
457 n=h_mesh%np; r=t; n=type; n=mode; n=n_start
subroutine, public init_velocity_pressure(mesh_f, mesh_c, time, dt, list_mode, un_m1, un, pn_m1, pn, phin_m1, phin)
real(kind=8) function, dimension(size(rr, 2)), public pp_exact(TYPE, rr, m, t)
subroutine error_petsc(string)
real(kind=8) function, dimension(size(rr, 2)), public vv_exact(TYPE, rr, m, t)
real(kind=8) function, dimension(size(rr, 2)), public temperature_exact(TYPE, rr, m, t)
real(kind=8) function, dimension(h_mesh%np), public extension_velocity(TYPE, H_mesh, mode, t, n_start)
real(kind=8) function, dimension(size(rr, 2)), public source_in_ns_momentum(TYPE, rr, mode, i, time, Re, ty, opt_density, opt_tempn)
subroutine, public init_temperature(mesh, time, dt, list_mode, tempn_m1, tempn)
real(kind=8) function, dimension(size(rr, 2)), public source_in_temperature(TYPE, rr, m, t)