Initialises the forcing data structures.
435 type(model_state_type),
target,
intent(inout) :: current_state
463 integer :: subsidence_input_type=divergence
468 character(len=STRING_LENGTH),
dimension(:),
allocatable :: units_q_force
469 character(len=STRING_LENGTH) :: units_theta_force=
'unset'
470 character(len=STRING_LENGTH) :: units_u_force=
'unset'
471 character(len=STRING_LENGTH) :: units_v_force=
'unset'
473 logical :: convert_input_theta_from_temperature=.false.
479 allocate(u_profile(current_state%local_grid%size(z_index)), &
480 v_profile(current_state%local_grid%size(z_index)), &
481 theta_profile(current_state%local_grid%size(z_index)), &
482 q_profile(current_state%local_grid%size(z_index)))
484 allocate(dtheta_profile(current_state%local_grid%size(z_index)), &
485 dq_profile(current_state%local_grid%size(z_index)), &
486 du_profile(current_state%local_grid%size(z_index)), &
487 dv_profile(current_state%local_grid%size(z_index)))
489 allocate(du_profile_diag(current_state%local_grid%size(z_index)), &
490 dv_profile_diag(current_state%local_grid%size(z_index)), &
491 dtheta_profile_diag(current_state%local_grid%size(z_index)), &
492 dq_profile_diag(current_state%local_grid%size(z_index), current_state%number_q_fields))
494 allocate(du_subs_profile_diag(current_state%local_grid%size(z_index)), &
495 dv_subs_profile_diag(current_state%local_grid%size(z_index)), &
496 dtheta_subs_profile_diag(current_state%local_grid%size(z_index)), &
497 dq_subs_profile_diag(current_state%local_grid%size(z_index), current_state%number_q_fields))
499 allocate(zgrid(current_state%local_grid%size(z_index)))
502 if (.not. current_state%passive_q .and. current_state%number_q_fields .gt. 0)
then
503 iqv=get_q_index(standard_q_names%VAPOUR,
'forcing')
504 iql=get_q_index(standard_q_names%CLOUD_LIQUID_MASS,
'forcing')
506 if (current_state%rain_water_mixing_ratio_index > 0) &
507 iqr = current_state%rain_water_mixing_ratio_index
508 if (current_state%ice_water_mixing_ratio_index > 0) &
509 iqi = current_state%ice_water_mixing_ratio_index
510 if (current_state%snow_water_mixing_ratio_index > 0) &
511 iqs = current_state%snow_water_mixing_ratio_index
512 if (current_state%graupel_water_mixing_ratio_index > 0) &
513 iqg = current_state%graupel_water_mixing_ratio_index
517 use_time_varying_subsidence= &
518 options_get_logical(current_state%options_database,
"use_time_varying_subsidence")
519 l_subs_pl_theta=options_get_logical(current_state%options_database,
"l_subs_pl_theta")
520 l_subs_pl_q=options_get_logical(current_state%options_database,
"l_subs_pl_q")
521 subsidence_input_type=options_get_integer(current_state%options_database,
"subsidence_input_type")
522 l_subs_local_theta=options_get_logical(current_state%options_database,
"subsidence_local_theta")
523 l_subs_local_q=options_get_logical(current_state%options_database,
"subsidence_local_q")
525 input_file=options_get_string(current_state%options_database,
"forcing_file")
527 if ((l_subs_pl_theta .and. .not. l_subs_local_theta) .or. &
528 (l_subs_pl_q .and. .not. l_subs_local_q))
then
529 if (.not. is_component_enabled(current_state%options_database,
"mean_profiles"))
then
534 if (l_subs_pl_theta .or. l_subs_pl_q)
then
536 if (trim(input_file)==
'' .or. trim(input_file)==
'None')
then
537 if (.not. use_time_varying_subsidence)
then
538 allocate(z_subs_pl(options_get_array_size(current_state%options_database,
"z_subs_pl")), &
539 f_subs_pl(options_get_array_size(current_state%options_database,
"f_subs_pl")))
540 call options_get_real_array(current_state%options_database,
"z_subs_pl", z_subs_pl)
541 call options_get_real_array(current_state%options_database,
"f_subs_pl", f_subs_pl)
543 zgrid=current_state%global_grid%configuration%vertical%z(:)
544 call piecewise_linear_1d(z_subs_pl(1:
size(z_subs_pl)), f_subs_pl(1:
size(f_subs_pl)), zgrid, &
545 current_state%global_grid%configuration%vertical%w_subs)
546 if (subsidence_input_type==divergence)
then
547 current_state%global_grid%configuration%vertical%w_subs(:) = &
548 -1.0*current_state%global_grid%configuration%vertical%w_subs(:)*zgrid(:)
553 deallocate(z_subs_pl, f_subs_pl)
555 if (use_time_varying_subsidence)
then
556 call check_forcing_status(nf90_open(path = trim(input_file), mode = nf90_nowrite, ncid = ncid))
561 call read_2d_forcing_dimensions(ncid, time_dim,z_dim)
562 allocate(forcing_input_times(time_dim), z_subs_pl(z_dim), f_subs_2d(z_dim, time_dim))
563 call read_2d_forcing_variables(trim(input_file), ncid, time_dim, forcing_input_times, &
564 z_dim, z_subs_pl, wsubs_key, f_subs_2d)
565 call check_forcing_status(nf90_close(ncid))
567 zgrid=current_state%global_grid%configuration%vertical%z(:)
569 allocate(current_state%global_grid%configuration%vertical%wsubs_time_vary(
size(zgrid), time_dim))
570 call piecewise_linear_2d(z_subs_pl(1:z_dim), forcing_input_times(1:time_dim), &
571 f_subs_2d(1:z_dim,1:time_dim), zgrid, &
572 current_state%global_grid%configuration%vertical%wsubs_time_vary)
581 if (.not.
allocated(current_state%l_forceq))
then
582 allocate(current_state%l_forceq(current_state%number_q_fields))
583 current_state%l_forceq=.false.
586 l_constant_forcing_theta=options_get_logical(current_state%options_database,
"l_constant_forcing_theta")
587 l_constant_forcing_q=options_get_logical(current_state%options_database,
"l_constant_forcing_q")
588 l_constant_forcing_u=options_get_logical(current_state%options_database,
"l_constant_forcing_u")
589 l_constant_forcing_v=options_get_logical(current_state%options_database,
"l_constant_forcing_v")
591 if (l_constant_forcing_q)
then
592 allocate(names_force_pl_q(options_get_array_size(current_state%options_database,
"names_constant_forcing_q")))
593 call options_get_string_array(current_state%options_database,
"names_constant_forcing_q", names_force_pl_q)
596 if (l_constant_forcing_theta)
then
597 constant_forcing_type_theta=options_get_integer(current_state%options_database,
"constant_forcing_type_theta")
598 forcing_timescale_theta=options_get_real(current_state%options_database,
"forcing_timescale_theta")
599 l_constant_forcing_theta_z2pressure=options_get_logical(current_state%options_database,
"l_constant_forcing_theta_z2pressure")
601 allocate(z_force_pl_theta(options_get_array_size(current_state%options_database,
"z_force_pl_theta")), &
602 f_force_pl_theta(options_get_array_size(current_state%options_database,
"f_force_pl_theta")))
603 call options_get_real_array(current_state%options_database,
"z_force_pl_theta", z_force_pl_theta)
604 call options_get_real_array(current_state%options_database,
"f_force_pl_theta", f_force_pl_theta)
606 relax_to_initial_theta_profile=options_get_logical(current_state%options_database,
"relax_to_initial_theta_profile")
607 if (relax_to_initial_theta_profile)
then
608 current_state%global_grid%configuration%vertical%theta_force(:) = &
609 current_state%global_grid%configuration%vertical%theta_init(:)
611 if (l_constant_forcing_theta_z2pressure)
then
612 zgrid=current_state%global_grid%configuration%vertical%zn(:)
614 zgrid=current_state%global_grid%configuration%vertical%prefn(:)
616 call piecewise_linear_1d(z_force_pl_theta(1:
size(z_force_pl_theta)), f_force_pl_theta(1:
size(f_force_pl_theta)), zgrid, &
617 current_state%global_grid%configuration%vertical%theta_force)
621 convert_input_theta_from_temperature=options_get_logical(current_state%options_database, &
622 "convert_input_theta_from_temperature")
623 if (convert_input_theta_from_temperature)
then
624 current_state%global_grid%configuration%vertical%theta_force(:) = &
625 current_state%global_grid%configuration%vertical%theta_force(:)* &
626 current_state%global_grid%configuration%vertical%prefrcp(:)
629 if (constant_forcing_type_theta==tendency)
then
630 units_theta_force=options_get_string(current_state%options_database,
"units_theta_force")
631 select case(trim(units_theta_force))
633 current_state%global_grid%configuration%vertical%theta_force(:) = &
634 current_state%global_grid%configuration%vertical%theta_force(:)/seconds_in_a_day
638 deallocate(z_force_pl_theta, f_force_pl_theta)
642 if (l_constant_forcing_u)
then
643 constant_forcing_type_u=options_get_integer(current_state%options_database,
"constant_forcing_type_u")
644 forcing_timescale_u=options_get_real(current_state%options_database,
"forcing_timescale_u")
645 relax_to_initial_u_profile=options_get_logical(current_state%options_database,
"relax_to_initial_u_profile")
646 if (relax_to_initial_u_profile)
then
647 current_state%global_grid%configuration%vertical%u_force(:) = &
648 current_state%global_grid%configuration%vertical%u_init(:)
650 allocate(z_force_pl_u(options_get_array_size(current_state%options_database,
"z_force_pl_u")), &
651 f_force_pl_u(options_get_array_size(current_state%options_database,
"f_force_pl_u")))
652 call options_get_real_array(current_state%options_database,
"z_force_pl_u", z_force_pl_u)
653 call options_get_real_array(current_state%options_database,
"f_force_pl_u", f_force_pl_u)
655 zgrid=current_state%global_grid%configuration%vertical%zn(:)
656 call piecewise_linear_1d(z_force_pl_u(1:
size(z_force_pl_u)), f_force_pl_u(1:
size(f_force_pl_u)), zgrid, &
657 current_state%global_grid%configuration%vertical%u_force)
658 deallocate(z_force_pl_u, f_force_pl_u)
662 if (constant_forcing_type_u==tendency)
then
664 units_u_force=options_get_string(current_state%options_database,
"units_u_force")
665 select case(trim(units_u_force))
666 case(m_per_second_per_day)
667 current_state%global_grid%configuration%vertical%u_force(:) = &
668 current_state%global_grid%configuration%vertical%u_force(:)/seconds_in_a_day
676 if (l_constant_forcing_v)
then
677 constant_forcing_type_v=options_get_integer(current_state%options_database,
"constant_forcing_type_v")
678 forcing_timescale_v=options_get_real(current_state%options_database,
"forcing_timescale_v")
679 relax_to_initial_v_profile=options_get_logical(current_state%options_database,
"relax_to_initial_v_profile")
680 if (relax_to_initial_v_profile)
then
681 current_state%global_grid%configuration%vertical%v_force(:) = &
682 current_state%global_grid%configuration%vertical%v_init(:)
684 allocate(z_force_pl_v(options_get_array_size(current_state%options_database,
"z_force_pl_v")), &
685 f_force_pl_v(options_get_array_size(current_state%options_database,
"f_force_pl_v")))
686 call options_get_real_array(current_state%options_database,
"z_force_pl_v", z_force_pl_v)
687 call options_get_real_array(current_state%options_database,
"f_force_pl_v", f_force_pl_v)
689 zgrid=current_state%global_grid%configuration%vertical%zn(:)
690 call piecewise_linear_1d(z_force_pl_v(1:
size(z_force_pl_v)), f_force_pl_v(1:
size(f_force_pl_v)), zgrid, &
691 current_state%global_grid%configuration%vertical%v_force)
692 deallocate(z_force_pl_v, f_force_pl_v)
696 if (constant_forcing_type_v==tendency)
then
698 units_v_force=options_get_string(current_state%options_database,
"units_v_force")
699 select case(trim(units_v_force))
700 case(m_per_second_per_day)
701 current_state%global_grid%configuration%vertical%v_force(:) = &
702 current_state%global_grid%configuration%vertical%v_force(:)/seconds_in_a_day
709 if (l_constant_forcing_q)
then
710 constant_forcing_type_q=options_get_integer(current_state%options_database,
"constant_forcing_type_q")
711 forcing_timescale_q=options_get_real(current_state%options_database,
"forcing_timescale_q")
712 nq_force=
size(names_force_pl_q)
713 allocate(z_force_pl_q(options_get_array_size(current_state%options_database,
"z_force_pl_q")))
714 call options_get_real_array(current_state%options_database,
"z_force_pl_q", z_force_pl_q)
715 nzq=
size(z_force_pl_q)
716 zgrid=current_state%global_grid%configuration%vertical%zn(:)
717 allocate(f_force_pl_q_tmp(nq_force*nzq))
718 call options_get_real_array(current_state%options_database,
"f_force_pl_q", f_force_pl_q_tmp)
719 allocate(f_force_pl_q(nzq, nq_force))
720 f_force_pl_q(1:nzq, 1:nq_force)=reshape(f_force_pl_q_tmp, (/nzq, nq_force/))
722 allocate(units_q_force(options_get_array_size(current_state%options_database,
"units_q_force")))
723 call options_get_string_array(current_state%options_database,
"units_q_force", units_q_force)
725 iq=get_q_index(trim(names_force_pl_q(n)),
'forcing:time-independent')
726 call piecewise_linear_1d(z_force_pl_q(1:nzq), f_force_pl_q(1:nzq,n), zgrid, &
727 current_state%global_grid%configuration%vertical%q_force(:,iq))
729 current_state%l_forceq(iq)=.true.
732 if (constant_forcing_type_q==tendency)
then
733 select case(trim(units_q_force(n)))
734 case(kg_per_kg_per_day)
735 current_state%global_grid%configuration%vertical%q_force(:,iq) = &
736 current_state%global_grid%configuration%vertical%q_force(:,iq)/seconds_in_a_day
737 case(g_per_kg_per_day)
738 current_state%global_grid%configuration%vertical%q_force(:,iq) = &
739 0.001*current_state%global_grid%configuration%vertical%q_force(:,iq)/seconds_in_a_day
740 case(g_per_kg_per_second)
741 current_state%global_grid%configuration%vertical%q_force(:,iq) = &
742 0.001*current_state%global_grid%configuration%vertical%q_force(:,iq)
747 deallocate(f_force_pl_q_tmp, units_q_force, f_force_pl_q, z_force_pl_q)
755 l_qdiag = (.not. current_state%passive_q .and. current_state%number_q_fields .gt. 0)
757 l_tend_pr_tot_u = current_state%u%active
758 l_tend_pr_tot_v = current_state%v%active
759 l_tend_pr_tot_th = current_state%th%active
760 l_tend_pr_tot_qv = l_qdiag .and. current_state%number_q_fields .ge. 1
761 l_tend_pr_tot_ql = l_qdiag .and. current_state%number_q_fields .ge. 2
762 l_tend_pr_tot_qi = l_qdiag .and. current_state%number_q_fields .ge. 11
763 l_tend_pr_tot_qr = l_qdiag .and. current_state%number_q_fields .ge. 11
764 l_tend_pr_tot_qs = l_qdiag .and. current_state%number_q_fields .ge. 11
765 l_tend_pr_tot_qg = l_qdiag .and. current_state%number_q_fields .ge. 11
766 l_tend_pr_tot_tabs = l_tend_pr_tot_th
768 l_tend_3d_u = current_state%u%active .or. l_tend_pr_tot_u
769 l_tend_3d_v = current_state%v%active .or. l_tend_pr_tot_v
770 l_tend_3d_th = current_state%th%active .or. l_tend_pr_tot_th
771 l_tend_3d_qv = (l_qdiag .and. current_state%number_q_fields .ge. 1) .or. l_tend_pr_tot_qv
772 l_tend_3d_ql = (l_qdiag .and. current_state%number_q_fields .ge. 2) .or. l_tend_pr_tot_ql
773 l_tend_3d_qi = (l_qdiag .and. current_state%number_q_fields .ge. 11) .or. l_tend_pr_tot_qi
774 l_tend_3d_qr = (l_qdiag .and. current_state%number_q_fields .ge. 11) .or. l_tend_pr_tot_qr
775 l_tend_3d_qs = (l_qdiag .and. current_state%number_q_fields .ge. 11) .or. l_tend_pr_tot_qs
776 l_tend_3d_qg = (l_qdiag .and. current_state%number_q_fields .ge. 11) .or. l_tend_pr_tot_qg
777 l_tend_3d_tabs = l_tend_3d_th
780 if (l_tend_3d_u)
then
781 allocate( tend_3d_u(current_state%local_grid%size(z_index), &
782 current_state%local_grid%size(y_index), &
783 current_state%local_grid%size(x_index) ) )
785 if (l_tend_3d_v)
then
786 allocate( tend_3d_v(current_state%local_grid%size(z_index), &
787 current_state%local_grid%size(y_index), &
788 current_state%local_grid%size(x_index) ) )
790 if (l_tend_3d_th)
then
791 allocate( tend_3d_th(current_state%local_grid%size(z_index), &
792 current_state%local_grid%size(y_index), &
793 current_state%local_grid%size(x_index) ) )
795 if (l_tend_3d_qv)
then
796 allocate( tend_3d_qv(current_state%local_grid%size(z_index), &
797 current_state%local_grid%size(y_index), &
798 current_state%local_grid%size(x_index) ) )
800 if (l_tend_3d_ql)
then
801 allocate( tend_3d_ql(current_state%local_grid%size(z_index), &
802 current_state%local_grid%size(y_index), &
803 current_state%local_grid%size(x_index) ) )
805 if (l_tend_3d_qi)
then
806 allocate( tend_3d_qi(current_state%local_grid%size(z_index), &
807 current_state%local_grid%size(y_index), &
808 current_state%local_grid%size(x_index) ) )
810 if (l_tend_3d_qr)
then
811 allocate( tend_3d_qr(current_state%local_grid%size(z_index), &
812 current_state%local_grid%size(y_index), &
813 current_state%local_grid%size(x_index) ) )
815 if (l_tend_3d_qs)
then
816 allocate( tend_3d_qs(current_state%local_grid%size(z_index), &
817 current_state%local_grid%size(y_index), &
818 current_state%local_grid%size(x_index) ) )
820 if (l_tend_3d_qg)
then
821 allocate( tend_3d_qg(current_state%local_grid%size(z_index), &
822 current_state%local_grid%size(y_index), &
823 current_state%local_grid%size(x_index) ) )
825 if (l_tend_3d_tabs)
then
826 allocate( tend_3d_tabs(current_state%local_grid%size(z_index), &
827 current_state%local_grid%size(y_index), &
828 current_state%local_grid%size(x_index) ) )
832 if (l_tend_pr_tot_u)
then
833 allocate( tend_pr_tot_u(current_state%local_grid%size(z_index)) )
835 if (l_tend_pr_tot_v)
then
836 allocate( tend_pr_tot_v(current_state%local_grid%size(z_index)) )
838 if (l_tend_pr_tot_th)
then
839 allocate( tend_pr_tot_th(current_state%local_grid%size(z_index)) )
841 if (l_tend_pr_tot_qv)
then
842 allocate( tend_pr_tot_qv(current_state%local_grid%size(z_index)) )
844 if (l_tend_pr_tot_ql)
then
845 allocate( tend_pr_tot_ql(current_state%local_grid%size(z_index)) )
847 if (l_tend_pr_tot_qi)
then
848 allocate( tend_pr_tot_qi(current_state%local_grid%size(z_index)) )
850 if (l_tend_pr_tot_qr)
then
851 allocate( tend_pr_tot_qr(current_state%local_grid%size(z_index)) )
853 if (l_tend_pr_tot_qs)
then
854 allocate( tend_pr_tot_qs(current_state%local_grid%size(z_index)) )
856 if (l_tend_pr_tot_qg)
then
857 allocate( tend_pr_tot_qg(current_state%local_grid%size(z_index)) )
859 if (l_tend_pr_tot_tabs)
then
860 allocate( tend_pr_tot_tabs(current_state%local_grid%size(z_index)) )
864 diagnostic_generation_frequency=options_get_integer(current_state%options_database,
"sampling_frequency")