73 type(model_state_type),
target,
intent(inout) :: current_state
74 character(len=*),
intent(in) :: name
75 type(component_field_information_type),
intent(out) :: field_information
79 field_information%field_type=component_array_field_type
80 field_information%data_type=component_double_data_type
81 field_information%number_dimensions=1
82 field_information%dimension_sizes(1)=current_state%local_grid%size(z_index)
83 field_information%enabled=.true.
86 strcomp=index(name,
"viscosity_3d_local")
87 if (strcomp .ne. 0)
then
88 field_information%field_type=component_array_field_type
89 field_information%number_dimensions=3
90 field_information%dimension_sizes(1)=current_state%local_grid%size(z_index)
91 field_information%dimension_sizes(2)=current_state%local_grid%size(y_index)
92 field_information%dimension_sizes(3)=current_state%local_grid%size(x_index)
93 field_information%data_type=component_double_data_type
95 if (name .eq.
"tend_u_viscosity_3d_local")
then
97 else if (name .eq.
"tend_v_viscosity_3d_local")
then
99 else if (name .eq.
"tend_w_viscosity_3d_local")
then
102 field_information%enabled=.true.
108 strcomp=index(name,
"viscosity_profile_total_local")
109 if (strcomp .ne. 0)
then
110 field_information%field_type=component_array_field_type
111 field_information%number_dimensions=1
112 field_information%dimension_sizes(1)=current_state%local_grid%size(z_index)
113 field_information%data_type=component_double_data_type
115 if (name .eq.
"tend_u_viscosity_profile_total_local")
then
117 else if (name .eq.
"tend_v_viscosity_profile_total_local")
then
119 else if (name .eq.
"tend_w_viscosity_profile_total_local")
then
122 field_information%enabled=.true.
135 type(model_state_type),
target,
intent(inout) :: current_state
136 character(len=*),
intent(in) :: name
137 type(component_field_value_type),
intent(out) :: field_value
139 if (name .eq.
"u_viscosity")
then
141 else if (name .eq.
"v_viscosity")
then
143 else if (name .eq.
"w_viscosity")
then
147 else if (name .eq.
"tend_u_viscosity_3d_local" .and.
allocated(
tend_3d_u))
then
149 else if (name .eq.
"tend_v_viscosity_3d_local" .and.
allocated(
tend_3d_v))
then
151 else if (name .eq.
"tend_w_viscosity_3d_local" .and.
allocated(
tend_3d_w))
then
155 else if (name .eq.
"tend_u_viscosity_profile_total_local" .and.
allocated(
tend_pr_tot_u))
then
157 else if (name .eq.
"tend_v_viscosity_profile_total_local" .and.
allocated(
tend_pr_tot_v))
then
159 else if (name .eq.
"tend_w_viscosity_profile_total_local" .and.
allocated(
tend_pr_tot_w))
then
169 type(model_state_type),
target,
intent(inout) :: current_state
171 integer :: z_size, y_size, x_size
173 z_size=current_state%local_grid%size(z_index) + current_state%local_grid%halo_size(z_index) * 2
174 y_size=current_state%local_grid%size(y_index) + current_state%local_grid%halo_size(y_index) * 2
175 x_size=current_state%local_grid%size(x_index) + current_state%local_grid%halo_size(x_index) * 2
176 allocate(current_state%vis_coefficient%data(z_size, y_size, x_size))
178 z_size=current_state%global_grid%size(z_index)
192 allocate(
tend_3d_u(current_state%local_grid%size(z_index), &
193 current_state%local_grid%size(y_index), &
194 current_state%local_grid%size(x_index) ) )
197 allocate(
tend_3d_v(current_state%local_grid%size(z_index), &
198 current_state%local_grid%size(y_index), &
199 current_state%local_grid%size(x_index) ) )
202 allocate(
tend_3d_w(current_state%local_grid%size(z_index), &
203 current_state%local_grid%size(y_index), &
204 current_state%local_grid%size(x_index) ) )
209 allocate(
tend_pr_tot_u(current_state%local_grid%size(z_index)) )
212 allocate(
tend_pr_tot_v(current_state%local_grid%size(z_index)) )
215 allocate(
tend_pr_tot_w(current_state%local_grid%size(z_index)) )
224 type(model_state_type),
target,
intent(inout) :: current_state
243 type(model_state_type),
target,
intent(inout) :: current_state
245 integer :: local_y, locaL_x, k, target_x_index, target_y_index
246 real(kind=default_precision),
dimension(current_state%local_grid%size(Z_INDEX)) :: tau12, tau12_ym1, tau12m1, &
247 tau11, tau22, tau22_yp1, tau33, tau23_ym1, tau11p1, tau13, tau13m1, tau23
249 local_y=current_state%column_local_y
250 local_x=current_state%column_local_x
251 target_y_index=local_y-current_state%local_grid%halo_size(y_index)
252 target_x_index=local_x-current_state%local_grid%halo_size(x_index)
255 if (current_state%first_timestep_column)
then
267 if (.not. current_state%use_viscosity_and_diffusion .or. current_state%halo_column)
return
269 if (current_state%viscosity_halo_swap_state%swap_in_progress)
then
271 call complete_nonblocking_halo_swap(current_state, current_state%viscosity_halo_swap_state, &
279 if (current_state%field_stepping == forward_stepping)
then
280 call calculate_tau(current_state, local_y, local_x, current_state%u, current_state%v, current_state%w, tau12, tau12_ym1, &
281 tau12m1, tau11, tau22, tau22_yp1, tau33, tau11p1, tau13, tau13m1, tau23, tau23_ym1)
283 call calculate_tau(current_state, local_y, local_x, current_state%zu, current_state%zv, current_state%zw, tau12, tau12_ym1, &
284 tau12m1, tau11, tau22, tau22_yp1, tau33, tau11p1, tau13, tau13m1, tau23, tau23_ym1)
286 call calculate_viscous_sources(current_state, local_y, local_x, tau12, tau12_ym1, tau12m1, tau11, tau22, tau22_yp1, tau33, &
287 tau11p1, tau13, tau13m1, tau23, tau23_ym1)
300 tau11, tau22, tau22_yp1, tau33, tau11p1, tau13, tau13m1, tau23, tau23_ym1)
301 type(model_state_type),
target,
intent(inout) :: current_state
302 integer,
intent(in) :: local_y, local_x
303 real(kind=default_precision),
dimension(:),
intent(in) :: tau12, tau12_ym1, tau12m1, tau11, tau22, tau22_yp1, tau33, &
304 tau11p1, tau13, tau13m1, tau23, tau23_ym1
308 do k=2, current_state%local_grid%size(z_index)
310 u_viscosity(k)=((tau11p1(k)-tau11(k))*current_state%global_grid%configuration%horizontal%cx+(tau12(k)-tau12_ym1(k))*&
311 current_state%global_grid%configuration%horizontal%cy+(tau13(k)-tau13(k-1))*&
312 current_state%global_grid%configuration%vertical%rdz(k))/current_state%global_grid%configuration%vertical%rhon(k)
313 current_state%su%data(k, local_y, local_x)=current_state%su%data(k, local_y, local_x)+
u_viscosity(k)
316 v_viscosity(k)=((tau12(k)-tau12m1(k))*current_state%global_grid%configuration%horizontal%cx+(tau22_yp1(k)-tau22(k))*&
317 current_state%global_grid%configuration%horizontal%cy+(tau23(k)-tau23(k-1))*&
318 current_state%global_grid%configuration%vertical%rdz(k))/current_state%global_grid%configuration%vertical%rhon(k)
319 current_state%sv%data(k, local_y, local_x)=current_state%sv%data(k, local_y, local_x)+
v_viscosity(k)
323 do k=2, current_state%local_grid%size(z_index)-1
324 w_viscosity(k)=((tau13(k)-tau13m1(k))*current_state%global_grid%configuration%horizontal%cx+(tau23(k)-tau23_ym1(k))*&
325 current_state%global_grid%configuration%horizontal%cy+(tau33(k+1)-tau33(k))*&
326 current_state%global_grid%configuration%vertical%rdzn(k+1))/current_state%global_grid%configuration%vertical%rho(k)
327 current_state%sw%data(k, local_y, local_x)=current_state%sw%data(k, local_y, local_x)+
w_viscosity(k)
336 subroutine calculate_tau(current_state, local_y, local_x, zu, zv, zw, tau12, tau12_ym1, tau12m1, tau11, &
337 tau22, tau22_yp1, tau33, tau11p1, tau13, tau13m1, tau23, tau23_ym1)
338 type(model_state_type),
target,
intent(inout) :: current_state
339 type(prognostic_field_type),
intent(inout) :: zu, zv, zw
340 integer,
intent(in) :: local_y, local_x
341 real(kind=default_precision),
dimension(:),
intent(out) :: tau12, tau12m1, tau11, tau22, tau22_yp1, tau33, &
342 tau11p1, tau13, tau13m1, tau23, tau23_ym1, tau12_ym1
345 real(kind=default_precision) :: vistmp, vis12, vis12a, visonp2, visonp2a, vis13, vis23
348 do k=1, current_state%local_grid%size(z_index)
350 vistmp=current_state%vis_coefficient%data(k, local_y, local_x)+current_state%vis_coefficient%data(k, local_y+1, local_x)+&
351 current_state%vis_coefficient%data(k-1, local_y, local_x)+current_state%vis_coefficient%data(k-1, local_y+1, local_x)
352 vis12=0.125_default_precision*(vistmp + current_state%vis_coefficient%data(k, local_y, local_x+1)+&
353 current_state%vis_coefficient%data(k, local_y+1, local_x+1)+&
354 current_state%vis_coefficient%data(k-1, local_y, local_x+1)+&
355 current_state%vis_coefficient%data(k-1, local_y+1, local_x+1))
356 tau12(k)=0.0_default_precision
358 tau12(k)=(zu%data(k, local_y+1, local_x)-zu%data(k, local_y, local_x))*&
359 current_state%global_grid%configuration%horizontal%cy
362 tau12(k)=tau12(k)+(zv%data(k, local_y, local_x+1)-zv%data(k, local_y, local_x))*&
363 current_state%global_grid%configuration%horizontal%cx
365 tau12(k)=tau12(k)*current_state%global_grid%configuration%vertical%rhon(k)*vis12
367 vis12a=0.125_default_precision*(vistmp +current_state%vis_coefficient%data(k, local_y, local_x-1)+&
368 current_state%vis_coefficient%data(k, local_y+1, local_x-1)+&
369 current_state%vis_coefficient%data(k-1, local_y, local_x-1)+&
370 current_state%vis_coefficient%data(k-1, local_y+1, local_x-1))
371 tau12m1(k)=0.0_default_precision
373 tau12m1(k)=(zu%data(k, local_y+1, local_x-1)-zu%data(k, local_y, local_x-1))*&
374 current_state%global_grid%configuration%horizontal%cy
377 tau12m1(k)=tau12m1(k)+(zv%data(k, local_y, local_x)-zv%data(k, local_y, local_x-1))*&
378 current_state%global_grid%configuration%horizontal%cx
380 tau12m1(k)=tau12m1(k)*current_state%global_grid%configuration%vertical%rhon(k)*vis12a
382 vistmp=current_state%vis_coefficient%data(k, local_y-1, local_x)+current_state%vis_coefficient%data(k, local_y, local_x)+&
383 current_state%vis_coefficient%data(k-1, local_y-1, local_x)+current_state%vis_coefficient%data(k-1, local_y, local_x)
384 vis12=0.125_default_precision*(vistmp + current_state%vis_coefficient%data(k, local_y-1, local_x+1)+&
385 current_state%vis_coefficient%data(k, local_y, local_x+1)+&
386 current_state%vis_coefficient%data(k-1, local_y-1, local_x+1)+&
387 current_state%vis_coefficient%data(k-1, local_y, local_x+1))
388 tau12_ym1(k)=0.0_default_precision
390 tau12_ym1(k)=(zu%data(k, local_y, local_x)-zu%data(k, local_y-1, local_x))*&
391 current_state%global_grid%configuration%horizontal%cy
394 tau12_ym1(k)=tau12_ym1(k)+(zv%data(k, local_y-1, local_x+1)-zv%data(k, local_y-1, local_x))*&
395 current_state%global_grid%configuration%horizontal%cx
397 tau12_ym1(k)=tau12_ym1(k)*current_state%global_grid%configuration%vertical%rhon(k)*vis12
399 visonp2=current_state%global_grid%configuration%vertical%rhon(k)*&
400 (current_state%vis_coefficient%data(k, local_y, local_x)+current_state%vis_coefficient%data(k-1, local_y, local_x))
402 tau11(k)=visonp2*(zu%data(k, local_y, local_x)-zu%data(k, local_y, local_x-1))*&
403 current_state%global_grid%configuration%horizontal%cx
405 tau11(k)=0.0_default_precision
408 tau22(k)=visonp2*(zv%data(k, local_y, local_x)-zv%data(k, local_y-1, local_x))*&
409 current_state%global_grid%configuration%horizontal%cy
411 tau22(k)=0.0_default_precision
414 tau33(k)=visonp2*(zw%data(k, local_y, local_x)-zw%data(k-1, local_y, local_x))*&
415 current_state%global_grid%configuration%vertical%rdz(k)
417 tau33(k)=0.0_default_precision
421 visonp2=current_state%global_grid%configuration%vertical%rhon(k)*&
422 (current_state%vis_coefficient%data(k, local_y+1, local_x)+&
423 current_state%vis_coefficient%data(k-1, local_y+1, local_x))
424 tau22_yp1(k)=visonp2*(zv%data(k, local_y+1, local_x)-zv%data(k, local_y, local_x))*&
425 current_state%global_grid%configuration%horizontal%cy
429 visonp2a=current_state%global_grid%configuration%vertical%rhon(k)*&
430 (current_state%vis_coefficient%data(k, local_y, local_x+1)+&
431 current_state%vis_coefficient%data(k-1, local_y, local_x+1))
432 tau11p1(k)=visonp2a*(current_state%zu%data(k, local_y, local_x+1)-&
433 zu%data(k, local_y, local_x))*current_state%global_grid%configuration%horizontal%cx
436 tau12(k)=0.0_default_precision
437 tau12_ym1(k)=0.0_default_precision
438 tau12m1(k)=0.0_default_precision
439 tau11(k)=0.0_default_precision
440 tau22(k)=0.0_default_precision
441 tau22_yp1(k)=0.0_default_precision
442 tau33(k)=0.0_default_precision
443 tau11p1(k)=0.0_default_precision
445 if (k .lt. current_state%local_grid%size(z_index))
then
446 vis13=0.5_default_precision*(current_state%vis_coefficient%data(k, local_y, local_x)+&
447 current_state%vis_coefficient%data(k, local_y, local_x+1))
448 tau13(k)=0.0_default_precision
450 tau13(k)=(zu%data(k+1, local_y, local_x)-zu%data(k, local_y, local_x))*&
451 current_state%global_grid%configuration%vertical%rdzn(k+1)
454 tau13(k)=tau13(k)+(zw%data(k, local_y, local_x+1)-zw%data(k, local_y, local_x))*&
455 current_state%global_grid%configuration%horizontal%cx
457 tau13(k)=tau13(k)*current_state%global_grid%configuration%vertical%rho(k)*vis13
459 vis13=0.5_default_precision*(current_state%vis_coefficient%data(k, local_y, local_x-1)+&
460 current_state%vis_coefficient%data(k, local_y, local_x))
461 tau13m1(k)=0.0_default_precision
463 tau13m1(k)=(zu%data(k+1, local_y, local_x-1)-zu%data(k, local_y, local_x-1))*&
464 current_state%global_grid%configuration%vertical%rdzn(k+1)
467 tau13m1(k)=tau13m1(k)+(zw%data(k, local_y, local_x)-zw%data(k, local_y, local_x-1))*&
468 current_state%global_grid%configuration%horizontal%cx
470 tau13m1(k)=tau13m1(k)*current_state%global_grid%configuration%vertical%rho(k)*vis13
472 vis23=0.5_default_precision*(current_state%vis_coefficient%data(k, local_y, local_x)+&
473 current_state%vis_coefficient%data(k, local_y+1, local_x))
474 tau23(k)=0.0_default_precision
476 tau23(k)=(zw%data(k, local_y+1, local_x)-zw%data(k, local_y, local_x))*&
477 current_state%global_grid%configuration%horizontal%cy
480 tau23(k)=tau23(k)+(zv%data(k+1, local_y, local_x)-zv%data(k, local_y, local_x))*&
481 current_state%global_grid%configuration%vertical%rdzn(k+1)
483 tau23(k)=tau23(k)*current_state%global_grid%configuration%vertical%rho(k)*vis23
485 vis23=0.5_default_precision*(current_state%vis_coefficient%data(k, local_y-1, local_x)+&
486 current_state%vis_coefficient%data(k, local_y, local_x))
487 tau23_ym1(k)=0.0_default_precision
489 tau23_ym1(k)=(zw%data(k, local_y, local_x)-zw%data(k, local_y-1, local_x))*&
490 current_state%global_grid%configuration%horizontal%cy
493 tau23_ym1(k)=tau23_ym1(k)+(zv%data(k+1, local_y-1, local_x)-zv%data(k, local_y-1, local_x))*&
494 current_state%global_grid%configuration%vertical%rdzn(k+1)
496 tau23_ym1(k)=tau23_ym1(k)*current_state%global_grid%configuration%vertical%rho(k)*vis23
498 tau13(k)=0.0_default_precision
499 tau13m1(k)=0.0_default_precision
500 tau23(k)=0.0_default_precision
501 tau23_ym1(k)=0.0_default_precision
510 type(model_state_type),
intent(inout) :: current_state
511 integer,
intent(in) :: halo_depth
512 logical,
intent(in) :: involve_corners
513 type(field_data_wrapper_type),
dimension(:),
intent(in),
optional :: source_data
515 call perform_local_data_copy_for_field(current_state%vis_coefficient%data, current_state%local_grid, &
516 current_state%parallel%my_rank, halo_depth, involve_corners)
528 neighbour_location, current_page, source_data)
529 type(model_state_type),
intent(inout) :: current_state
530 integer,
intent(in) :: dim, target_index, neighbour_location
531 integer,
intent(inout) :: current_page(:)
532 type(neighbour_description_type),
intent(inout) :: neighbour_description
533 type(field_data_wrapper_type),
dimension(:),
intent(in),
optional :: source_data
535 call copy_buffer_to_field(current_state%local_grid, neighbour_description%recv_halo_buffer, &
536 current_state%vis_coefficient%data, dim, target_index, current_page(neighbour_location))
538 current_page(neighbour_location)=current_page(neighbour_location)+1
551 y_target_index, neighbour_location, current_page, source_data)
552 type(model_state_type),
intent(inout) :: current_state
553 integer,
intent(in) :: corner_loc, x_target_index, y_target_index, neighbour_location
554 integer,
intent(inout) :: current_page(:)
555 type(neighbour_description_type),
intent(inout) :: neighbour_description
556 type(field_data_wrapper_type),
dimension(:),
intent(in),
optional :: source_data
558 call copy_buffer_to_corner(current_state%local_grid, neighbour_description%recv_corner_buffer, &
559 current_state%vis_coefficient%data, corner_loc, x_target_index, y_target_index, current_page(neighbour_location))
561 current_page(neighbour_location)=current_page(neighbour_location)+1
571 type(model_state_type),
target,
intent(in) :: current_state
572 integer,
intent(in) :: cxn, cyn, txn, tyn
576 tend_3d_u(:,tyn,txn)=current_state%su%data(:,cyn,cxn)
579 tend_3d_v(:,tyn,txn)=current_state%sv%data(:,cyn,cxn)
582 tend_3d_w(:,tyn,txn)=current_state%sw%data(:,cyn,cxn)
595 type(model_state_type),
target,
intent(inout) :: current_state
596 integer,
intent(in) :: cxn, cyn, txn, tyn
628 type(component_field_value_type),
intent(inout) :: field_value
629 real(kind=default_precision),
dimension(:),
optional :: real_1d_field
630 real(kind=default_precision),
dimension(:,:),
optional :: real_2d_field
631 real(kind=default_precision),
dimension(:,:,:),
optional :: real_3d_field
633 if (
present(real_1d_field))
then
634 allocate(field_value%real_1d_array(
size(real_1d_field)), source=real_1d_field)
635 else if (
present(real_2d_field))
then
636 allocate(field_value%real_2d_array(
size(real_2d_field, 1),
size(real_2d_field, 2)), source=real_2d_field)
637 else if (
present(real_3d_field))
then
638 allocate(field_value%real_3d_array(
size(real_3d_field, 1),
size(real_3d_field, 2),
size(real_3d_field, 3)), &
639 source=real_3d_field)