diff --git a/CONTRIBUTORS.md b/CONTRIBUTORS.md index 15cc3962b..97c5a022d 100644 --- a/CONTRIBUTORS.md +++ b/CONTRIBUTORS.md @@ -10,4 +10,5 @@ | harry-shepherd | Harry Shepherd | Met Office | 2026-01-08 | | DrTVockerodtMO | Terence Vockerodt | Met Office | 2026-01-08 | | MetBenjaminWent | Benjamin Went | Met Office | 2026-01-15 | -| timgraham-Met | Tim Graham | Met Office | 2026-01-15 | \ No newline at end of file +| timgraham-Met | Tim Graham | Met Office | 2026-01-15 | +| cjohnson-pi | Christine Johnson | Met Office | 2026-01-16 | \ No newline at end of file diff --git a/applications/adjoint_tests/example/configuration.nml b/applications/adjoint_tests/example/configuration.nml index 8e7dd58f1..645cbd2f3 100644 --- a/applications/adjoint_tests/example/configuration.nml +++ b/applications/adjoint_tests/example/configuration.nml @@ -180,8 +180,10 @@ write_fluxes=.false., write_minmax_tseries=.false., / &linear +fixed_ls=.false. pert_option='analytic', l_stabilise_bl=.false., +transport_efficiency=.false. / &logging run_log_level='info', @@ -297,11 +299,11 @@ panel_edge_treatment='none' profile_size=5, reversible=5*.false., runge_kutta_method='ssp3', -scheme=5*1, +scheme=5*3, si_outer_transport='none', slice_order='parabola', special_edges_monotone=0,0,0,1,1 -splitting=5*1, +splitting=5*2, substep_transport='off', theta_dispersion_correction=.false., theta_variable='dry', diff --git a/applications/adjoint_tests/source/algorithm/transport/common/adjt_end_transport_step_alg_mod.x90 b/applications/adjoint_tests/source/algorithm/transport/common/adjt_end_transport_step_alg_mod.x90 index 9342ce993..7fa6bf71e 100644 --- a/applications/adjoint_tests/source/algorithm/transport/common/adjt_end_transport_step_alg_mod.x90 +++ b/applications/adjoint_tests/source/algorithm/transport/common/adjt_end_transport_step_alg_mod.x90 @@ -21,7 +21,7 @@ module adjt_end_transport_step_alg_mod use model_clock_mod, only : model_clock_type use finite_element_config_mod, only : element_order_h, & element_order_v - use fs_continuity_mod, only : W2 + use fs_continuity_mod, only : W2, W2V use function_space_collection_mod, only : function_space_collection use setup_test_alg_mod, only : setup_test_tl_transport_controller use init_from_controller_alg_mod, only : init_counter_fieldvals @@ -78,6 +78,7 @@ module adjt_end_transport_step_alg_mod ! Variables for initialising fields type(function_space_type), pointer :: vector_space_w2_ptr + type(function_space_type), pointer :: vector_space_w2v_ptr ! Inner products real(kind=r_tran) :: inner1 @@ -108,13 +109,16 @@ module adjt_end_transport_step_alg_mod vector_space_w2_ptr => function_space_collection%get_fs( & mesh, element_order_h, element_order_v, W2 & - ) + ) + vector_space_w2v_ptr => function_space_collection%get_fs( & + mesh, element_order_h, element_order_v, W2V) + ! Assume that the last split step is in W2V (for VHV splitting) call sum_flux%initialise( vector_space = vector_space_w2_ptr ) - call flux_last_step%initialise( vector_space = vector_space_w2_ptr ) + call flux_last_step%initialise( vector_space = vector_space_w2v_ptr ) call sum_flux_input%initialise( vector_space = vector_space_w2_ptr ) - call flux_last_step_input%initialise( vector_space = vector_space_w2_ptr ) + call flux_last_step_input%initialise( vector_space = vector_space_w2v_ptr ) sum_flux_inner_prod = 0.0_r_def flux_last_step_inner_prod = 0.0_r_def diff --git a/applications/adjoint_tests/source/algorithm/transport/common/atlt_end_transport_step_alg_mod.x90 b/applications/adjoint_tests/source/algorithm/transport/common/atlt_end_transport_step_alg_mod.x90 index 08a3b7a2b..a4d64e064 100644 --- a/applications/adjoint_tests/source/algorithm/transport/common/atlt_end_transport_step_alg_mod.x90 +++ b/applications/adjoint_tests/source/algorithm/transport/common/atlt_end_transport_step_alg_mod.x90 @@ -21,12 +21,14 @@ module atlt_end_transport_step_alg_mod use model_clock_mod, only : model_clock_type use finite_element_config_mod, only : element_order_h, & element_order_v - use fs_continuity_mod, only : W2, W3 + use fs_continuity_mod, only : W2, W3, W2V use function_space_collection_mod, only : function_space_collection use tl_transport_controller_mod, only : tl_transport_controller_type use transport_controller_mod, only : transport_controller_type use flux_precomputations_alg_mod, only : flux_precomputations_type use transport_counter_mod, only : transport_counter_type + use split_transport_utils_mod, only : get_num_split_steps, & + get_direction_w2_fs use setup_test_alg_mod, only : setup_test_tl_transport_controller use transport_metadata_collection_mod, only : transport_metadata_collection use transport_metadata_mod, only : transport_metadata_type @@ -233,6 +235,9 @@ module atlt_end_transport_step_alg_mod type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_mass_2_inp ! Variables for initialising fields + integer(kind=i_def) :: splitting + integer(kind=i_def) :: num_split_steps + integer(kind=i_def) :: fs_id type(function_space_type), pointer :: vector_space_w2_ptr type(function_space_type), pointer :: vector_space_w3_ptr @@ -284,11 +289,16 @@ module atlt_end_transport_step_alg_mod dummy_ref_mass_2_inp, & pert_transport_controller ) - vector_space_w2_ptr => function_space_collection%get_fs( & - mesh, element_order_h, element_order_v, W2 & + ! W2 function space appropriate for splitting (fs on the last step) + splitting = transport_metadata%get_splitting() + num_split_steps = get_num_split_steps(splitting) + fs_id = get_direction_w2_fs(splitting, num_split_steps) + vector_space_w2_ptr => function_space_collection%get_fs( & + mesh, element_order_h, element_order_v, fs_id & ) - vector_space_w3_ptr => function_space_collection%get_fs( & - mesh, element_order_h, element_order_v, W3 & + + vector_space_w3_ptr => function_space_collection%get_fs( & + mesh, element_order_h, element_order_v, W3 & ) call field_np1%initialise( vector_space = vector_space_w3_ptr ) diff --git a/applications/adjoint_tests/source/algorithm/transport/control/atlt_transport_field_alg_mod.x90 b/applications/adjoint_tests/source/algorithm/transport/control/atlt_transport_field_alg_mod.x90 index d1120743d..6ec2d229b 100644 --- a/applications/adjoint_tests/source/algorithm/transport/control/atlt_transport_field_alg_mod.x90 +++ b/applications/adjoint_tests/source/algorithm/transport/control/atlt_transport_field_alg_mod.x90 @@ -27,13 +27,19 @@ module atlt_transport_field_alg_mod use setup_test_alg_mod, only : setup_test_tl_transport_controller use inner_from_controller_rdef_alg_mod, only : flux_pc_ref_flux_prod_rdef, & flux_pc_ref_flux_inp_prod_rdef, & + counter_flux_prod_rdef, & + counter_flux_inp_prod_rdef, & + counter_field_n_prod_rdef, & + counter_field_n_inp_prod_rdef, & wind_pc_dir_prod_rdef, & wind_pc_dir_inp_prod_rdef use init_from_controller_alg_mod, only : init_flux_pc_fieldvals, & - init_wind_pc_fieldvals + init_wind_pc_fieldvals, & + init_counter_fieldvals use transport_enumerated_types_mod, only : direction_3d, & direction_h, & direction_v + use fs_continuity_mod, only : W3 implicit none @@ -86,18 +92,25 @@ module atlt_transport_field_alg_mod ! Variables used to handle calculations for fields stored in tl_transport_controller integer(kind=i_def) :: mesh_id + type(transport_controller_type), pointer :: transport_controller type(transport_controller_type), pointer :: pert_transport_controller type(transport_controller_type), pointer :: ls_transport_controller type(transport_metadata_type), pointer :: transport_metadata type(r_tran_field_type), dimension(:), allocatable :: fpc_ls_wind_ref_flux_inp + type(r_tran_field_type), dimension(:), allocatable :: fpc_ls_wind_ref_flux_inp2 type(r_tran_field_type), dimension(:), allocatable :: fpc_pert_wind_ref_flux_inp type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_field_1_inp type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_mass_1_inp type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_field_2_inp type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_mass_2_inp + type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_field_3_inp + type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_mass_3_inp type(r_tran_field_type) :: wind_pc_h_inp type(r_tran_field_type) :: wind_pc_v_inp type(r_tran_field_type) :: wind_pc_3d_inp + type(r_tran_field_type) :: tctr_field_n_inp + type(r_tran_field_type), dimension(:), allocatable :: tctr_flux_inp + type(function_space_type), pointer :: vector_space_w3_ptr ! Test parameters and variables real(kind=r_def), parameter :: overall_tolerance = 1500.0_r_def @@ -106,11 +119,16 @@ module atlt_transport_field_alg_mod mesh_id = mesh%get_id() call setup_test_tl_transport_controller( mesh, model_clock, tl_transport_controller, field_n, field_n_fs ) + transport_controller => tl_transport_controller%get_ls_wind_ls_rho_controller() ls_transport_controller => tl_transport_controller%get_ls_wind_pert_rho_controller() pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() transport_metadata => pert_transport_controller%get_transport_metadata() ! Initialising fields + call init_counter_fieldvals( mesh, & + tctr_field_n_inp, & + tctr_flux_inp, & + pert_transport_controller ) call init_flux_pc_fieldvals( mesh, & 1_i_def, & fpc_ls_wind_ref_flux_inp, & @@ -123,6 +141,12 @@ module atlt_transport_field_alg_mod dummy_ref_field_2_inp, & dummy_ref_mass_2_inp, & pert_transport_controller ) + call init_flux_pc_fieldvals( mesh, & + 1_i_def, & + fpc_ls_wind_ref_flux_inp2, & + dummy_ref_field_3_inp, & + dummy_ref_mass_3_inp, & + transport_controller ) call init_wind_pc_fieldvals( mesh, & wind_pc_h_inp, & direction_h, & @@ -136,6 +160,11 @@ module atlt_transport_field_alg_mod direction_3d, & pert_transport_controller ) + vector_space_w3_ptr => function_space_collection%get_fs( & + mesh, element_order_h, element_order_v, W3 & + ) + call field_n%initialise(vector_space = vector_space_w3_ptr) + call field_n%copy_field_properties( field_np1 ) call field_n%copy_field_properties( ls_field_n ) @@ -147,6 +176,7 @@ module atlt_transport_field_alg_mod ! Initialise values and call the tangent-linear alg. call invoke( setval_random( field_np1 ), & + setval_random( field_n), & setval_x( field_np1_input, field_np1 ), & setval_x( field_n_input, field_n ), & setval_random( ls_field_n ) ) @@ -179,6 +209,7 @@ module atlt_transport_field_alg_mod call wind_pc_dir_prod_rdef( inner1, pert_transport_controller, mesh_id, direction_h ) call wind_pc_dir_prod_rdef( inner1, pert_transport_controller, mesh_id, direction_v ) call wind_pc_dir_prod_rdef( inner1, pert_transport_controller, mesh_id, direction_3d ) + call counter_flux_prod_rdef( inner1, pert_transport_controller ) ! Scaling fields call invoke( inc_a_times_X( field_np1_sf, field_np1 ), & @@ -210,6 +241,7 @@ module atlt_transport_field_alg_mod call wind_pc_dir_inp_prod_rdef( inner2, pert_transport_controller, mesh_id, direction_h, wind_pc_h_inp ) call wind_pc_dir_inp_prod_rdef( inner2, pert_transport_controller, mesh_id, direction_v, wind_pc_v_inp ) call wind_pc_dir_inp_prod_rdef( inner2, pert_transport_controller, mesh_id, direction_3d, wind_pc_3d_inp ) + call counter_flux_inp_prod_rdef( inner2, pert_transport_controller, tctr_flux_inp ) call tl_transport_controller%finalise() if ( allocated(fpc_ls_wind_ref_flux_inp) ) deallocate( fpc_ls_wind_ref_flux_inp ) diff --git a/applications/adjoint_tests/source/algorithm/transport/init_from_controller_alg_mod.x90 b/applications/adjoint_tests/source/algorithm/transport/init_from_controller_alg_mod.x90 index 17e300fa3..d33548660 100644 --- a/applications/adjoint_tests/source/algorithm/transport/init_from_controller_alg_mod.x90 +++ b/applications/adjoint_tests/source/algorithm/transport/init_from_controller_alg_mod.x90 @@ -222,7 +222,9 @@ module init_from_controller_alg_mod if ( flux_ptr%is_initialised() ) then call invoke( setval_X( counter_flux_inp(step), flux_ptr ) ) else - call invoke( setval_c( counter_flux_inp(step), 0.0_r_tran ) ) + call counter_flux_inp(step)%copy_field_properties(flux_ptr) + call invoke( setval_random( flux_ptr ), & + setval_X( counter_flux_inp(step), flux_ptr ) ) end if end do diff --git a/applications/adjoint_tests/source/algorithm/transport/mol/atlt_mol_conservative_alg_mod.x90 b/applications/adjoint_tests/source/algorithm/transport/mol/atlt_mol_conservative_alg_mod.x90 index 6f2a2cabc..c33a61db0 100644 --- a/applications/adjoint_tests/source/algorithm/transport/mol/atlt_mol_conservative_alg_mod.x90 +++ b/applications/adjoint_tests/source/algorithm/transport/mol/atlt_mol_conservative_alg_mod.x90 @@ -97,16 +97,20 @@ module atlt_mol_conservative_alg_mod integer(kind=i_def) :: mesh_id type(transport_controller_type), pointer :: pert_transport_controller type(transport_controller_type), pointer :: ls_transport_controller + type(transport_controller_type), pointer :: transport_controller type(field_type) :: tctr_field_n type(function_space_type), pointer :: tctr_field_n_vs type(field_type) :: tctr_field_n_inp type(r_tran_field_type), dimension(:), allocatable :: tctr_flux_inp type(r_tran_field_type), dimension(:), allocatable :: fpc_ls_wind_ref_flux_inp type(r_tran_field_type), dimension(:), allocatable :: fpc_pert_wind_ref_flux_inp + type(r_tran_field_type), dimension(:), allocatable :: fpc_ls_wind_ref_flux_inp2 type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_field_1_inp type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_mass_1_inp type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_field_2_inp type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_mass_2_inp + type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_field_3_inp + type(r_tran_field_type), dimension(:), allocatable :: dummy_ref_mass_3_inp type(r_tran_field_type) :: wind_pc_h_inp type(r_tran_field_type) :: wind_pc_v_inp type(r_tran_field_type) :: wind_pc_3d_inp @@ -118,6 +122,7 @@ module atlt_mol_conservative_alg_mod mesh_id = mesh%get_id() call setup_test_tl_transport_controller( mesh, model_clock, tl_transport_controller, tctr_field_n, tctr_field_n_vs ) + transport_controller => tl_transport_controller%get_ls_wind_ls_rho_controller() ls_transport_controller => tl_transport_controller%get_ls_wind_pert_rho_controller() pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() @@ -126,6 +131,12 @@ module atlt_mol_conservative_alg_mod tctr_field_n_inp, & tctr_flux_inp, & pert_transport_controller ) + call init_flux_pc_fieldvals( mesh, & + 1_i_def, & + fpc_ls_wind_ref_flux_inp2, & + dummy_ref_field_3_inp, & + dummy_ref_mass_3_inp, & + transport_controller ) call init_flux_pc_fieldvals( mesh, & 1_i_def, & fpc_ls_wind_ref_flux_inp, & diff --git a/applications/linear_model/example/configuration.nml b/applications/linear_model/example/configuration.nml index 5ed21fc6c..f41da1efb 100644 --- a/applications/linear_model/example/configuration.nml +++ b/applications/linear_model/example/configuration.nml @@ -154,6 +154,7 @@ write_minmax_tseries=.false., fixed_ls=.false. pert_option='analytic', l_stabilise_bl=.false., +transport_efficiency=.false., / &logging run_log_level='info', @@ -262,11 +263,11 @@ panel_edge_high_order=.true., panel_edge_treatment='none' reversible=.false.,.false.,.false.,.false.,.false. runge_kutta_method='ssp3' -scheme=1,1,1,1,1 +scheme=5*3 si_outer_transport='none', slice_order='parabola' special_edges_monotone=5*1 -splitting=1,1,1,1,1 +splitting=5*2 substep_transport='off' theta_dispersion_correction=.false. theta_variable='dry' diff --git a/rose-stem/app/adjoint_tests/rose-app.conf b/rose-stem/app/adjoint_tests/rose-app.conf index 20c70408b..dfe65d41b 100644 --- a/rose-stem/app/adjoint_tests/rose-app.conf +++ b/rose-stem/app/adjoint_tests/rose-app.conf @@ -1333,11 +1333,11 @@ panel_edge_treatment='none' profile_size=5 reversible=5*.false. runge_kutta_method='ssp3' -scheme=5*1 +scheme=5*3 si_outer_transport='none' slice_order='parabola' special_edges_monotone=5*1 -splitting=5*1 +splitting=5*2 substep_transport='off' theta_dispersion_correction=.false. theta_variable='dry' diff --git a/rose-stem/app/linear_model/opt/rose-app-runge-kutta.conf b/rose-stem/app/linear_model/opt/rose-app-runge-kutta.conf index 1047582b2..4f831da45 100644 --- a/rose-stem/app/linear_model/opt/rose-app-runge-kutta.conf +++ b/rose-stem/app/linear_model/opt/rose-app-runge-kutta.conf @@ -98,4 +98,6 @@ cfl_mol_2d_stab=2.0 cfl_mol_3d_stab=2.0 max_vert_cfl_calc='uniform' runge_kutta_method='ssp4' +scheme=5*1 +splitting=5*1 !!wind_mono_top_depth=0 diff --git a/rose-stem/app/linear_model/opt/rose-app-semi-implicit.conf b/rose-stem/app/linear_model/opt/rose-app-semi-implicit.conf index be23c9c84..b300e53d1 100644 --- a/rose-stem/app/linear_model/opt/rose-app-semi-implicit.conf +++ b/rose-stem/app/linear_model/opt/rose-app-semi-implicit.conf @@ -87,3 +87,5 @@ tau_u=0.5 [namelist:transport] !!wind_mono_top_depth=0 +scheme=5*1 +splitting=5*1 \ No newline at end of file diff --git a/rose-stem/app/linear_model/rose-app.conf b/rose-stem/app/linear_model/rose-app.conf index 1c83c6635..78d13e797 100644 --- a/rose-stem/app/linear_model/rose-app.conf +++ b/rose-stem/app/linear_model/rose-app.conf @@ -725,6 +725,7 @@ ls_read_w2h=.false. max_bl_stabilisation=0.75 n_bl_levels_to_stabilise=15 pert_option='file' +transport_efficiency=.true. [namelist:logging] log_to_rank_zero_only=.false. @@ -1376,11 +1377,11 @@ panel_edge_treatment='none' profile_size=5 reversible=.true.,.true.,.false.,.true.,.true. runge_kutta_method='ssp3' -scheme=5*1 +scheme=5*3 si_outer_transport='none' slice_order='parabola' special_edges_monotone=5*1 -splitting=5*1 +splitting=5*2 substep_transport='off' theta_dispersion_correction=.false. theta_variable='dry' diff --git a/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_nwp_gal9-C12_azspice_gnu_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_nwp_gal9-C12_azspice_gnu_fast-debug-64bit.txt index b17a42b79..007c1751c 100644 --- a/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_nwp_gal9-C12_azspice_gnu_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_nwp_gal9-C12_azspice_gnu_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 3FA7A1197E0CB9A0 -Inner product checksum theta = 4161FFC69BF205C3 -Inner product checksum u = 4569B9484E575B5E -Inner product checksum mr1 = 3F36C0A882D88E66 -Inner product checksum mr2 = 3EEADFD85794EA13 -Inner product checksum mr3 = 3EDFEC0C984D5704 -Inner product checksum mr4 = 3EED3A950512EB1A +Inner product checksum rho = 3FA79A9D0614EEEE +Inner product checksum theta = 4162215D8ED8A125 +Inner product checksum u = 456A083CF999EC55 +Inner product checksum mr1 = 3F36D11B3C65CF3F +Inner product checksum mr2 = 3EEAE32B46C09D7B +Inner product checksum mr3 = 3EDFEC26C0DD3244 +Inner product checksum mr4 = 3EED3B0E06CBD627 Inner product checksum mr5 = 0 Inner product checksum mr6 = 0 diff --git a/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_azspice_gnu_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_azspice_gnu_fast-debug-64bit.txt index 36256062e..7d6d62d97 100644 --- a/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_azspice_gnu_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_azspice_gnu_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 3F95A1DDD2EC7839 -Inner product checksum theta = 414000E886CBFB32 -Inner product checksum u = 452CCCFAA0BBF41E -Inner product checksum mr1 = 3F2D2183BE8E5C80 -Inner product checksum mr2 = 3EDF53A45356ADEA -Inner product checksum mr3 = 3ED26AEF2D8453B5 -Inner product checksum mr4 = 3EE009842E269AF6 +Inner product checksum rho = 3F95ACC94441A964 +Inner product checksum theta = 413FD47788AE2B0E +Inner product checksum u = 452CBAADFF659374 +Inner product checksum mr1 = 3F2D499F1DCEDE0B +Inner product checksum mr2 = 3EDF578E9B65451B +Inner product checksum mr3 = 3ED26AFC441BB781 +Inner product checksum mr4 = 3EE00986BD98837B Inner product checksum mr5 = 0 Inner product checksum mr6 = 0 diff --git a/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_op_azspice_gnu_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_op_azspice_gnu_fast-debug-64bit.txt index 5ae16e5aa..fc254a564 100644 --- a/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_op_azspice_gnu_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/jedi_lfric_tests/azspice/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_op_azspice_gnu_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 3F95A1DF1E431E10 -Inner product checksum theta = 41400115CD2072D9 -Inner product checksum u = 452CCC1C9C56D3B6 -Inner product checksum mr1 = 3F2D21B9057A6428 -Inner product checksum mr2 = 3EDF53C45279F22E -Inner product checksum mr3 = 3ED26AEE8FE0B97C -Inner product checksum mr4 = 3EE009839EC7E10F +Inner product checksum rho = 3F95ACCAC59D64D9 +Inner product checksum theta = 413FD499BD5603B7 +Inner product checksum u = 452CB9ABAE30ABAF +Inner product checksum mr1 = 3F2D49D33C196D4C +Inner product checksum mr2 = 3EDF57AEC293C509 +Inner product checksum mr3 = 3ED26AFB9145F0F5 +Inner product checksum mr4 = 3EE009862BB0DCF3 Inner product checksum mr5 = 0 Inner product checksum mr6 = 0 diff --git a/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_nwp_gal9-C12_ex1a_cce_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_nwp_gal9-C12_ex1a_cce_fast-debug-64bit.txt index db90ee0d3..9ba16ef2b 100644 --- a/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_nwp_gal9-C12_ex1a_cce_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_nwp_gal9-C12_ex1a_cce_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 3FA7A119D15694EA -Inner product checksum theta = 4161FFC6B41D6C18 -Inner product checksum u = 4569B9484F679FDB -Inner product checksum mr1 = 3F36C0A85083444F -Inner product checksum mr2 = 3EEADFD859B930EE -Inner product checksum mr3 = 3EDFEC0C98392EF7 -Inner product checksum mr4 = 3EED3A95035B479C +Inner product checksum rho = 3FA79A9D24D80B84 +Inner product checksum theta = 4162215D998A790E +Inner product checksum u = 456A083CF563E471 +Inner product checksum mr1 = 3F36D11B19E12F6E +Inner product checksum mr2 = 3EEAE32B48FCBD56 +Inner product checksum mr3 = 3EDFEC26C094A93A +Inner product checksum mr4 = 3EED3B0E07D01F92 Inner product checksum mr5 = 0 Inner product checksum mr6 = 0 diff --git a/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_ex1a_cce_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_ex1a_cce_fast-debug-64bit.txt index 7dfc5de39..50ce4226d 100644 --- a/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_ex1a_cce_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_ex1a_cce_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 3F95A1C504B6A918 -Inner product checksum theta = 413FF942C2E563C0 -Inner product checksum u = 452CCCF6ABE0865F -Inner product checksum mr1 = 3F2D21C8E6687010 -Inner product checksum mr2 = 3EDF53AD2309BABA -Inner product checksum mr3 = 3ED26AEFEBB2B535 -Inner product checksum mr4 = 3EE00982584FBC51 +Inner product checksum rho = 3F95ACB014DA78D0 +Inner product checksum theta = 413FCBF88C4910CE +Inner product checksum u = 452CBA65B40AFAB8 +Inner product checksum mr1 = 3F2D49E974B2F080 +Inner product checksum mr2 = 3EDF5792D6D2D34E +Inner product checksum mr3 = 3ED26AFD01439F84 +Inner product checksum mr4 = 3EE00984EA466A8E Inner product checksum mr5 = 0 Inner product checksum mr6 = 0 diff --git a/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_op_ex1a_cce_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_op_ex1a_cce_fast-debug-64bit.txt index f6ca49a10..a8eadce8d 100644 --- a/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_op_ex1a_cce_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/jedi_lfric_tests/ex1a/checksum_jedi_lfric_tests_tlm_forecast_tl_default-C12_op_ex1a_cce_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 3F95A1E1F17F7A09 -Inner product checksum theta = 413FF80D6A5F4BE8 -Inner product checksum u = 452CCB73BB6342D4 -Inner product checksum mr1 = 3F2D21C0C423794A -Inner product checksum mr2 = 3EDF53B78E4C3283 -Inner product checksum mr3 = 3ED26AEF7B3222E8 -Inner product checksum mr4 = 3EE00982FE30ACFA +Inner product checksum rho = 3F95ACCD3EEAC2B8 +Inner product checksum theta = 413FCAAF700BA76A +Inner product checksum u = 452CB8D93BF67253 +Inner product checksum mr1 = 3F2D49E4C797F50E +Inner product checksum mr2 = 3EDF579D1DD1B948 +Inner product checksum mr3 = 3ED26AFC92438788 +Inner product checksum mr4 = 3EE00985A17972EA Inner product checksum mr5 = 0 Inner product checksum mr6 = 0 diff --git a/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_dcmip301-C24_azspice_gnu_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_dcmip301-C24_azspice_gnu_fast-debug-64bit.txt index 4148368a6..7e2da881a 100644 --- a/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_dcmip301-C24_azspice_gnu_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_dcmip301-C24_azspice_gnu_fast-debug-64bit.txt @@ -1,3 +1,3 @@ -Inner product checksum rho = 3F20250B4FE9B9EF -Inner product checksum theta = 403305E848D6E311 -Inner product checksum u = 433053B3FA026478 +Inner product checksum rho = 3F2022594816DB4A +Inner product checksum theta = 4033068A55A83CF8 +Inner product checksum u = 433056118E026D02 diff --git a/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_nwp_gal9-C12_MG_azspice_gnu_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_nwp_gal9-C12_MG_azspice_gnu_fast-debug-64bit.txt index 6a7be832b..421106c0a 100644 --- a/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_nwp_gal9-C12_MG_azspice_gnu_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_nwp_gal9-C12_MG_azspice_gnu_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 3FA958F0C1D55513 -Inner product checksum theta = 4161FFBE6B89106C -Inner product checksum u = 456A4D0AF00B4480 -Inner product checksum mr1 = 3F3B888F79EB2E18 -Inner product checksum mr2 = 3EEB3109C670AF3C -Inner product checksum mr3 = 3EDFF31B816930A0 -Inner product checksum mr4 = 3EED42827F73DFDC +Inner product checksum rho = 3FA958AB1C8F215C +Inner product checksum theta = 416354444EC30001 +Inner product checksum u = 456C296AB5699467 +Inner product checksum mr1 = 3F3BB706EC65A9C8 +Inner product checksum mr2 = 3EEB5F98BD833463 +Inner product checksum mr3 = 3EE049E45013AA3C +Inner product checksum mr4 = 3EED67E57E3C63B8 Inner product checksum mr5 = 0 Inner product checksum mr6 = 0 diff --git a/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_nwp_gal9_random-C12_MG_azspice_gnu_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_nwp_gal9_random-C12_MG_azspice_gnu_fast-debug-64bit.txt index c19cd8f5f..16ae3f01f 100644 --- a/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_nwp_gal9_random-C12_MG_azspice_gnu_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/linear_model/azspice/checksum_linear_model_nwp_gal9_random-C12_MG_azspice_gnu_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 40C2C1694FA30DD0 -Inner product checksum theta = 43076DE653AEC2FE -Inner product checksum u = 4697E3FAA368DE1C -Inner product checksum mr1 = 40A1D167F8ADA142 -Inner product checksum mr2 = 40699AF5BDBC4864 -Inner product checksum mr3 = 4069A3412D88D047 -Inner product checksum mr4 = 40697EB0A557078A +Inner product checksum rho = 40C231C005491FDA +Inner product checksum theta = 43076F9C99ACD554 +Inner product checksum u = 4697A0E8F7E333C4 +Inner product checksum mr1 = 40A1E7982CFF402D +Inner product checksum mr2 = 4069A4088DCC70F6 +Inner product checksum mr3 = 4069ABB4FA045ED2 +Inner product checksum mr4 = 406986D0FD3118C3 Inner product checksum mr5 = 406997D0C9B8EEFB Inner product checksum mr6 = 40699547929283AC diff --git a/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_dcmip301-C24_ex1a_gnu_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_dcmip301-C24_ex1a_gnu_fast-debug-64bit.txt index a5449d77a..489629753 100644 --- a/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_dcmip301-C24_ex1a_gnu_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_dcmip301-C24_ex1a_gnu_fast-debug-64bit.txt @@ -1,3 +1,3 @@ -Inner product checksum rho = 3F20250B300C3D30 -Inner product checksum theta = 403305E8470B6A22 -Inner product checksum u = 433053B3A41C6049 +Inner product checksum rho = 3F20225947FE2721 +Inner product checksum theta = 4033068A55825E07 +Inner product checksum u = 433056118DDB6C73 diff --git a/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_nwp_gal9-C12_MG_ex1a_gnu_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_nwp_gal9-C12_MG_ex1a_gnu_fast-debug-64bit.txt index 8839d63b9..6ad269013 100644 --- a/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_nwp_gal9-C12_MG_ex1a_gnu_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_nwp_gal9-C12_MG_ex1a_gnu_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 3FA958F09D97D52D -Inner product checksum theta = 4161FFBE7B6B7DE1 -Inner product checksum u = 456A4D0ACAB467CA -Inner product checksum mr1 = 3F3B888D5EA6A45C -Inner product checksum mr2 = 3EEB310977D9678E -Inner product checksum mr3 = 3EDFF31B705B97D8 -Inner product checksum mr4 = 3EED42824B13222B +Inner product checksum rho = 3FA958AC8E6AE778 +Inner product checksum theta = 4163544440909176 +Inner product checksum u = 456C296A75CB61E1 +Inner product checksum mr1 = 3F3BB70A4A2CDB2C +Inner product checksum mr2 = 3EEB5F9888ADF09C +Inner product checksum mr3 = 3EE049E45255EF3E +Inner product checksum mr4 = 3EED67E58C2021B8 Inner product checksum mr5 = 0 Inner product checksum mr6 = 0 diff --git a/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_nwp_gal9_random-C12_MG_ex1a_gnu_fast-debug-64bit.txt b/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_nwp_gal9_random-C12_MG_ex1a_gnu_fast-debug-64bit.txt index 1fba69899..a311bce1d 100644 --- a/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_nwp_gal9_random-C12_MG_ex1a_gnu_fast-debug-64bit.txt +++ b/rose-stem/site/meto/kgos/linear_model/ex1a/checksum_linear_model_nwp_gal9_random-C12_MG_ex1a_gnu_fast-debug-64bit.txt @@ -1,9 +1,9 @@ -Inner product checksum rho = 40C2C1668AA536C6 -Inner product checksum theta = 43076DE652800BE2 -Inner product checksum u = 4697E3FAA2E5E0C6 -Inner product checksum mr1 = 40A1D167F834C0E0 -Inner product checksum mr2 = 40699AF5BDC7A01E -Inner product checksum mr3 = 4069A3412D91ECF2 -Inner product checksum mr4 = 40697EB0A55E4C2B +Inner product checksum rho = 40C231C0161316BD +Inner product checksum theta = 43076F9C97F79AE5 +Inner product checksum u = 4697A0E8F3CD0B5C +Inner product checksum mr1 = 40A1E79829B25F82 +Inner product checksum mr2 = 4069A4088E1A7842 +Inner product checksum mr3 = 4069ABB4FA005079 +Inner product checksum mr4 = 406986D0FD2D0E7A Inner product checksum mr5 = 406997D0C9B8EEFC Inner product checksum mr6 = 40699547929283AD diff --git a/science/adjoint/source/algorithm/transport/common/adj_end_transport_step_alg_mod.x90 b/science/adjoint/source/algorithm/transport/common/adj_end_transport_step_alg_mod.x90 index eb4c4e7ba..3a35463d6 100644 --- a/science/adjoint/source/algorithm/transport/common/adj_end_transport_step_alg_mod.x90 +++ b/science/adjoint/source/algorithm/transport/common/adj_end_transport_step_alg_mod.x90 @@ -97,6 +97,14 @@ contains hori_flux => transport_counter%get_flux(2) call sum_vert_flux%initialise(w2v_fs) + if (.not. hori_flux%is_initialised() ) then + call hori_flux%initialise(w2h_fs) + call invoke( setval_c( hori_flux, 0.0_r_tran ) ) + end if + if (.not. old_flux%is_initialised() ) then + call old_flux%initialise(w2v_fs) + call invoke( setval_c( old_flux, 0.0_r_tran ) ) + end if call invoke( setval_c( sum_vert_flux, 0.0_r_tran ), & adj_combine_w2_field_kernel_type( sum_flux, hori_flux, & sum_vert_flux, & @@ -104,6 +112,7 @@ contains face_selector_ns ), & inc_X_plus_Y( old_flux, sum_vert_flux ), & inc_X_plus_Y( flux_last_step, sum_vert_flux ), & + setval_c( sum_flux, 0.0_r_tran ), & setval_c( sum_vert_flux, 0.0_r_tran ) ) ! General case for any splitting diff --git a/science/adjoint/source/algorithm/transport/common/adj_flux_precomputations_mod.x90 b/science/adjoint/source/algorithm/transport/common/adj_flux_precomputations_mod.x90 index f5c0dcb86..a742942de 100644 --- a/science/adjoint/source/algorithm/transport/common/adj_flux_precomputations_mod.x90 +++ b/science/adjoint/source/algorithm/transport/common/adj_flux_precomputations_mod.x90 @@ -141,7 +141,7 @@ module adj_flux_precomputations_mod class(flux_precomputations_type), intent(inout) :: flux_pc integer(kind=i_def), intent(in) :: step - type(r_tran_field_type), intent(in) :: ref_flux + type(r_tran_field_type), intent(inout) :: ref_flux ! Local variables integer(kind=i_def) :: splitting diff --git a/science/adjoint/source/algorithm/transport/common/atl_end_transport_step_alg_mod.x90 b/science/adjoint/source/algorithm/transport/common/atl_end_transport_step_alg_mod.x90 index 6f0cf3a45..ef1cfa971 100644 --- a/science/adjoint/source/algorithm/transport/common/atl_end_transport_step_alg_mod.x90 +++ b/science/adjoint/source/algorithm/transport/common/atl_end_transport_step_alg_mod.x90 @@ -38,7 +38,6 @@ module atl_end_transport_step_alg_mod use transport_metadata_mod, only: transport_metadata_type use adj_flux_precomputations_mod, only: adj_initialise_step - implicit none private @@ -223,14 +222,18 @@ contains if (.not. final_split_step) then flux_ptr => transport_counter%get_flux( step ) - call invoke( inc_X_plus_Y( flux, flux_ptr ), & - setval_c( flux_ptr, 0.0_r_tran ) ) + call invoke( inc_X_plus_Y( flux, flux_ptr ) ) end if ! Combine fluxes from each part of the calculation call invoke( inc_X_plus_Y( flux_pert_wind, flux ), & inc_X_plus_Y( flux_ls_wind, flux ) ) + if (.not. final_split_step) then + flux_ptr => transport_counter%get_flux( step ) + call invoke( setval_c( flux_ptr, 0.0_r_tran ) ) + end if + end subroutine atl_end_conservative_step_alg end module atl_end_transport_step_alg_mod diff --git a/science/adjoint/source/algorithm/transport/control/atl_split_transport_mod.x90 b/science/adjoint/source/algorithm/transport/control/atl_split_transport_mod.x90 new file mode 100644 index 000000000..2295ba4f1 --- /dev/null +++ b/science/adjoint/source/algorithm/transport/control/atl_split_transport_mod.x90 @@ -0,0 +1,334 @@ +!----------------------------------------------------------------------------- +! (C) Crown copyright 2026 Met Office. All rights reserved. +! The file LICENCE, distributed with this code, contains details of the terms +! under which the code may be used. +!----------------------------------------------------------------------------- +!> @brief Adjoint of the tangent linear version of split_transport + +module atl_split_transport_mod + + use constants_mod, only: i_def, r_tran + use r_tran_field_mod, only: r_tran_field_type + use log_mod, only: log_event, LOG_LEVEL_ERROR + use split_transport_utils_mod, only: get_splitting_direction, & + get_splitting_fraction, & + get_num_split_steps + use split_transport_mod, only: split_transport_field + + ! Transport control infrastructure + use transport_controller_mod, only: transport_controller_type + use tl_transport_controller_mod, only: tl_transport_controller_type + use transport_counter_mod, only: transport_counter_type + use transport_enumerated_types_mod, only: direction_h, & + direction_v, & + split_method_null, & + split_method_mol, & + split_method_ffsl, & + split_method_sl, & + equation_form_conservative, & + equation_form_advective, & + equation_form_consistent + use transport_metadata_mod, only: transport_metadata_type + + ! Algorithms + use end_of_transport_step_alg_mod, only: end_of_conservative_step_alg, & + end_of_consistent_step_alg, & + end_of_advective_step_alg + use linear_config_mod, only: transport_efficiency + + implicit none + + private + + public :: atl_split_transport_control + public :: atl_split_transport_field + +contains + + !============================================================================= + !> @brief Controls vertical/horizontal split transport. + !> @details Manages the vertical/horizontal splitting of the split transport + !! scheme by choosing the splitting type and calling the individual + !! vertical and horizontal split components. + !> @param[in,out] field_np1 ACTIVE Field to return at end + !> @param[in,out] field_n ACTIVE Field at the start + !> @param[in] ls_field_n PASSIVE LS Field at the start + !> @param[in,out] tl_transport_controller + !! Encapsulating object containing the + !! transport counter and precomputations + subroutine atl_split_transport_control(field_np1, field_n, ls_field_n, & + tl_transport_controller) + + implicit none + + ! Arguments + type(r_tran_field_type), intent(inout) :: field_np1 + type(r_tran_field_type), target, intent(inout) :: field_n + type(r_tran_field_type), target, intent(in) :: ls_field_n + type(tl_transport_controller_type), intent(inout) :: tl_transport_controller + + ! Internal variables + integer(kind=i_def) :: num_split_steps + integer(kind=i_def) :: split_step_count + type(r_tran_field_type), allocatable :: ls_field_np1(:) + type(r_tran_field_type), target :: field_tmp + type(r_tran_field_type), pointer :: field_ptr + type(transport_counter_type), pointer :: transport_counter + type(transport_counter_type), pointer :: ls_transport_counter + type(transport_metadata_type), pointer :: transport_metadata + type(transport_controller_type), pointer :: pert_transport_controller + type(transport_controller_type), pointer :: ls_transport_controller + real(kind=r_tran) :: dt_substep + integer(kind=i_def) :: num_substeps + + ! ------------------------------------------------------------------------ ! + ! NONLINEAR (LS) + ! ------------------------------------------------------------------------ ! + + ls_transport_controller => tl_transport_controller%get_ls_wind_ls_rho_controller() + transport_metadata => ls_transport_controller%get_transport_metadata() + ls_transport_counter => ls_transport_controller%get_transport_counter() + + ! Initialise fields + num_split_steps = get_num_split_steps(transport_metadata%get_splitting()) + allocate(ls_field_np1(num_split_steps + 1)) + do split_step_count = 1, num_split_steps + 1 + call field_n%copy_field_properties(ls_field_np1(split_step_count)) + end do + + if ( .not. transport_efficiency ) then + + ! When we have multiple split steps, we need an intermediate field for the + ! field at the start of each substep + ! field_ptr points to the field at the start of each split step + call invoke( setval_X(ls_field_np1(1), ls_field_n) ) + + if (num_split_steps > 1) then + call field_n%copy_field_properties(field_tmp) + call invoke( setval_X(field_tmp, ls_field_n) ) + field_ptr => field_tmp + else + field_ptr => ls_field_n + end if + + do split_step_count = 1, num_split_steps + call split_transport_field( ls_field_np1(split_step_count+1), & + field_ptr, & + ls_transport_controller & + ) + + if (split_step_count < num_split_steps) then + call invoke( setval_X(field_tmp, ls_field_np1(split_step_count+1) ) ) + + ! Increment split step counter + call ls_transport_counter%inc_split_step_counter() + end if + end do + else + do split_step_count = 1, num_split_steps + 1 + call invoke( setval_X(ls_field_np1(split_step_count), ls_field_n) ) + end do + end if + + ! ------------------------------------------------------------------------ ! + ! LINEAR + ! ------------------------------------------------------------------------ ! + pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() + transport_metadata => pert_transport_controller%get_transport_metadata() + transport_counter => pert_transport_controller%get_transport_counter() + + dt_substep = transport_counter%get_dt_substep() + num_substeps = transport_counter%get_num_substeps() + call transport_counter%adj_initialise(transport_metadata, dt_substep, num_substeps ) + call ls_transport_counter%adj_initialise(transport_metadata, dt_substep, num_substeps ) + + ! Initialise fields + num_split_steps = get_num_split_steps(transport_metadata%get_splitting()) + + call field_n%copy_field_properties(field_tmp) + + ! Set fields to zero + call invoke( setval_C(field_tmp, 0.0_r_def ) ) + + ! Loop backwards + do split_step_count = num_split_steps, 1, -1 + + if (split_step_count < num_split_steps) then + + ! Increment split step counter + call transport_counter%dec_split_step_counter() + call ls_transport_counter%dec_split_step_counter() + + call invoke( inc_X_plus_Y(field_np1, field_tmp ) ) + call invoke( setval_C(field_tmp, 0.0_r_def ) ) + + end if + + call atl_split_transport_field( & + field_np1, field_tmp, ls_field_np1(split_step_count), & + tl_transport_controller & + ) + call invoke( setval_C(field_np1, 0.0_r_def ) ) + + end do + + call invoke( setval_C(field_np1, 0.0_r_def ) ) + call invoke( inc_X_plus_Y(field_n, field_tmp) ) + call invoke( setval_C(field_tmp, 0.0_r_def ) ) + + deallocate( ls_field_np1 ) + + end subroutine atl_split_transport_control + + !============================================================================= + !> @brief Does either vertical or horizontal transport of a field. + !> @details Performs a vertical or horizontal transport step, solving the + !! transport equation for a (multidata) field. + !> @param[in,out] field_np1 ACTIVE Field to return at end + !> @param[in,out] field_n ACTIVE Field at the start + !> @param[in] ls_field_n PASSIVE Field at the start + !> @param[in,out] tl_transport_controller + !! Encapsulating object containing the + !! transport counter and precomputations + subroutine atl_split_transport_field(field_np1, field_n, ls_field_n, tl_transport_controller) + + use atl_mol_conservative_alg_mod, only: atl_mol_conservative_alg + use atl_mol_advective_alg_mod, only: atl_mol_advective_alg + + implicit none + + ! Arguments + type(r_tran_field_type), intent(inout) :: field_np1 + type(r_tran_field_type), intent(inout) :: field_n + type(r_tran_field_type), intent(in) :: ls_field_n + type(tl_transport_controller_type), intent(inout) :: tl_transport_controller + + ! Internal variables + integer(kind=i_def) :: method, direction + type(transport_counter_type), pointer :: transport_counter + type(transport_metadata_type), pointer :: transport_metadata + type(transport_controller_type), pointer :: pert_transport_controller + + pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() + transport_counter => pert_transport_controller%get_transport_counter() + transport_metadata => pert_transport_controller%get_transport_metadata() + + ! -------------------------------------------------------------------------! + ! Set up method based on direction + ! -------------------------------------------------------------------------! + direction = get_splitting_direction( & + transport_metadata%get_splitting(), & + transport_counter%get_split_step_of_substep_counter() & + ) + + select case ( direction ) + case ( direction_h ) + method = transport_metadata%get_horizontal_method() + case ( direction_v ) + method = transport_metadata%get_vertical_method() + case default + call log_event('Split transport direction not recognised', LOG_LEVEL_ERROR) + end select + + ! -------------------------------------------------------------------------! + ! Choose method, and then choose equation + ! -------------------------------------------------------------------------! + select case ( method ) + + ! -------------------------------------------------------------------------! + ! Null step + ! -------------------------------------------------------------------------! + case ( split_method_null ) + call log_event( & + 'ADJ: Split method null not coded', & + LOG_LEVEL_ERROR & + ) + + ! -------------------------------------------------------------------------! + ! Method of Lines step + ! -------------------------------------------------------------------------! + case ( split_method_mol ) + ! Choose form of transport equation + select case ( transport_metadata%get_equation_form() ) + case ( equation_form_conservative ) + call atl_mol_conservative_alg( & + field_np1, field_n, ls_field_n, tl_transport_controller & + ) + case ( equation_form_advective ) + call atl_mol_advective_alg( & + field_np1, field_n, ls_field_n, tl_transport_controller & + ) + case ( equation_form_consistent ) + call log_event( & + 'ADJ: MoL Consistent not coded yet', & + LOG_LEVEL_ERROR & + ) + case default + call log_event( & + 'Trying to solve unrecognised form of transport equation', & + LOG_LEVEL_ERROR & + ) + + end select + + ! -------------------------------------------------------------------------! + ! Flux-Form Semi-Lagrangian step + ! -------------------------------------------------------------------------! + case ( split_method_ffsl ) + ! All equation forms have the same control method + call log_event( & + 'ADJ: FFSL not coded yet', & + LOG_LEVEL_ERROR & + ) + ! -------------------------------------------------------------------------! + ! Semi-Lagrangian step + ! -------------------------------------------------------------------------! + case ( split_method_sl ) + ! Choose direction + select case ( direction ) + + case ( direction_h ) + ! Horizontal SL only for advective form + if ( transport_metadata%get_equation_form() /= equation_form_advective ) then + call log_event( & + 'Horizontal semi-Lagrangian is only for advective form', & + LOG_LEVEL_ERROR & + ) + end if + call log_event( & + 'ADJ: Horizontal SL not coded yet', & + LOG_LEVEL_ERROR & + ) + + case ( direction_v ) + ! Choose form of transport equation for vertical + select case ( transport_metadata%get_equation_form() ) + + case ( equation_form_conservative ) + call log_event( & + 'TL: Vertical SL conservative not coded yet', & + LOG_LEVEL_ERROR & + ) + case ( equation_form_advective ) + call log_event( & + 'ADJ: Vertical SL advective not coded yet', & + LOG_LEVEL_ERROR & + ) + case default + call log_event( & + 'Trying to solve unrecognised form of transport equation', & + LOG_LEVEL_ERROR & + ) + end select + end select + + case default + call log_event( & + 'Trying to transport with unrecognised scheme', & + LOG_LEVEL_ERROR & + ) + end select + + end subroutine atl_split_transport_field + +end module atl_split_transport_mod diff --git a/science/adjoint/source/algorithm/transport/control/atl_transport_field_mod.f90 b/science/adjoint/source/algorithm/transport/control/atl_transport_field_mod.f90 index a00c2ee0b..2a7a05cf9 100644 --- a/science/adjoint/source/algorithm/transport/control/atl_transport_field_mod.f90 +++ b/science/adjoint/source/algorithm/transport/control/atl_transport_field_mod.f90 @@ -20,6 +20,7 @@ module atl_transport_field_mod equation_form_advective use atl_mol_conservative_alg_mod, only: atl_mol_conservative_alg use atl_mol_advective_alg_mod, only: atl_mol_advective_alg + use atl_split_transport_mod, only: atl_split_transport_control use transport_controller_mod, only: transport_controller_type use tl_transport_controller_mod, only: tl_transport_controller_type use transport_counter_mod, only: transport_counter_type @@ -112,8 +113,8 @@ subroutine atl_transport_field(field_np1, field_n, ls_field_n, & ! Some split horizontal/vertical transport scheme ! -------------------------------------------------------------------------! case ( scheme_split ) - call log_event('Split transport not implemented for tangent-linear app', & - LOG_LEVEL_ERROR) + call atl_split_transport_control(field_np1, field_n, ls_field_n, & + tl_transport_controller) case default call log_event('Trying to transport with unrecognised scheme', & diff --git a/science/adjoint/source/algorithm/transport/mol/atl_mol_advective_alg_mod.x90 b/science/adjoint/source/algorithm/transport/mol/atl_mol_advective_alg_mod.x90 index cec84c386..3f63651f1 100644 --- a/science/adjoint/source/algorithm/transport/mol/atl_mol_advective_alg_mod.x90 +++ b/science/adjoint/source/algorithm/transport/mol/atl_mol_advective_alg_mod.x90 @@ -27,7 +27,8 @@ module atl_mol_advective_alg_mod use wind_precomputations_alg_mod, only: wind_precomputations_type ! Configuration - use transport_config_mod, only: runge_kutta_method + use transport_config_mod, only: runge_kutta_method + use linear_config_mod, only: transport_efficiency implicit none @@ -56,7 +57,9 @@ contains ! Internal variables integer(kind=i_def) :: stage, s, direction integer(kind=i_def) :: nstage, substep + integer(kind=i_def) :: ls_stage, ls_nstage integer(kind=i_def) :: number_substeps + integer(kind=i_def) :: ls_number_substeps, ls_substep integer(kind=i_def) :: step integer(kind=i_def) :: splitting real(kind=r_def) :: dt_mol_substep @@ -71,7 +74,10 @@ contains type(field_type), allocatable :: rk_field(:) real(kind=r_def), allocatable :: rk_weights(:,:) type(transport_metadata_type), pointer :: transport_metadata + type(transport_metadata_type), pointer :: ls_transport_metadata type(transport_counter_type), pointer :: transport_counter + type(transport_counter_type), pointer :: ls_transport_counter + type(transport_controller_type), pointer :: transport_controller type(transport_controller_type), pointer :: ls_transport_controller type(transport_controller_type), pointer :: pert_transport_controller type(wind_precomputations_type), pointer :: ls_wind_precomputations @@ -81,10 +87,16 @@ contains ! Extract transport objects and initialise temporary fields ! ------------------------------------------------------------------------ ! mesh => field%get_mesh() + transport_controller => tl_transport_controller%get_ls_wind_ls_rho_controller() ls_transport_controller => tl_transport_controller%get_ls_wind_pert_rho_controller() pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() + transport_metadata => pert_transport_controller%get_transport_metadata() + ls_transport_metadata => transport_controller%get_transport_metadata() + transport_counter => pert_transport_controller%get_transport_counter() + ls_transport_counter => transport_controller%get_transport_counter() + ls_wind_precomputations => ls_transport_controller%get_wind_precomputations() pert_wind_precomputations => pert_transport_controller%get_wind_precomputations() @@ -102,13 +114,24 @@ contains ) call get_rk_transport_weights(nstage, rk_weights, runge_kutta_method) + + if (transport_efficiency) then + ls_number_substeps = 1 + ls_nstage =1 + else + ls_number_substeps = number_substeps + ls_nstage = nstage + end if + allocate( rk_field(nstage) ) - allocate( stored_ls_field(number_substeps, nstage) ) + allocate( stored_ls_field(ls_number_substeps, ls_nstage) ) do stage = 1, nstage call rk_field(stage)%initialise(field%get_function_space()) - do substep = 1, number_substeps - call stored_ls_field(substep, stage)%initialise(field%get_function_space()) + end do + do ls_stage = 1, ls_nstage + do ls_substep = 1, ls_number_substeps + call stored_ls_field(ls_substep, ls_stage)%initialise(field%get_function_space()) end do end do call rhs_field%initialise( field%get_function_space() ) @@ -135,19 +158,23 @@ contains ! array for each substep and each RK stage to use in the perturbation !--------------------------------------------------------------------------! + if (transport_efficiency) then + call invoke( setval_X(stored_ls_field(1,1), ls_field_np1 )) + else + ! Perform the number of rk-stages and substeps required - do substep = 1, number_substeps + do substep = 1, ls_number_substeps ! Reset field_n ready for the this substep call ls_field_np1%copy_field_properties(ls_field_n) call invoke( setval_X(ls_field_n, ls_field_np1) ) - do stage = 1, nstage + do stage = 1, ls_nstage ! Store values for use in the perturbation code. This is done before ! advecting so as to pick up the correct linearisation state call invoke( setval_X(stored_ls_field(substep, stage), ls_field_np1) ) - final_rk_stage = ( stage == nstage ) + final_rk_stage = ( stage == ls_nstage ) ! Compute the field for this stage: ! rhs_field = sum(s=1,stage): a(stage,s)*field^(s) @@ -158,9 +185,9 @@ contains end do ! Compute update: rhs = u.grad(rhs_field) - call advective_and_flux_alg(dummy, rhs, rhs_field, ls_field_n, & - ls_advecting_wind, direction, & - transport_metadata, final_rk_stage, & + call advective_and_flux_alg(dummy, rhs, rhs_field, ls_field_n, & + ls_advecting_wind, direction, & + ls_transport_metadata, final_rk_stage, & dt_mol_substep, .false., .true. ) ! Update field: f = f^n - dt*rhs @@ -170,10 +197,12 @@ contains end do ! End of step: if necessary enforce min val and overwrite in blending zone - call end_of_advective_step_alg( & - ls_field_np1, ls_field_n, transport_counter, transport_metadata & + call end_of_advective_step_alg( & + ls_field_np1, ls_field_n, ls_transport_counter, ls_transport_metadata & ) + end if + ! -------------------------------------------------------------------------- ! Perturbation ! @@ -199,8 +228,20 @@ contains ! Perform the number of rk-stages and substeps required do substep = number_substeps, 1, -1 + if (transport_efficiency) then + ls_substep = 1 + else + ls_substep = substep + end if + do stage = nstage, 1, -1 + if (transport_efficiency) then + ls_stage = 1 + else + ls_stage = stage + end if + final_rk_stage = ( stage == nstage ) ! Update field: f = f^n - dt_substep*rhs @@ -211,7 +252,7 @@ contains ! Compute update: rhs = u.grad(rhs_field) call atl_advective_and_flux_alg( & dummy, dummy, rhs, rhs_field, advecting_wind, & - stored_ls_field(substep, stage), ls_advecting_wind, & + stored_ls_field(ls_substep, ls_stage), ls_advecting_wind, & direction, transport_metadata, final_rk_stage, dt_mol_substep, & .false., .true. & ) diff --git a/science/adjoint/source/algorithm/transport/mol/atl_mol_conservative_alg_mod.x90 b/science/adjoint/source/algorithm/transport/mol/atl_mol_conservative_alg_mod.x90 index b5ed18d7e..75bd9d794 100644 --- a/science/adjoint/source/algorithm/transport/mol/atl_mol_conservative_alg_mod.x90 +++ b/science/adjoint/source/algorithm/transport/mol/atl_mol_conservative_alg_mod.x90 @@ -26,7 +26,7 @@ module atl_mol_conservative_alg_mod use runge_kutta_init_mod, only: get_rk_transport_weights use split_transport_utils_mod, only: get_num_split_steps, & get_splitting_direction - use transport_constants_mod, only: get_directional_im3_div + use transport_constants_mod, only: get_directional_im3_div_r_tran use transport_controller_mod, only: transport_controller_type use tl_transport_controller_mod, only: tl_transport_controller_type use transport_counter_mod, only: transport_counter_type @@ -40,6 +40,7 @@ module atl_mol_conservative_alg_mod ! Configuration use boundaries_config_mod, only: limited_area use base_mesh_config_mod, only: topology, topology_non_periodic + use linear_config_mod, only: transport_efficiency use io_config_mod, only: subroutine_timers use transport_config_mod, only: runge_kutta_method, & dry_field_name, & @@ -74,7 +75,9 @@ module atl_mol_conservative_alg_mod integer(kind=i_def) :: mesh_id integer(kind=i_def) :: stage, s, direction integer(kind=i_def) :: nstage, substep - integer(kind=i_def) :: number_substeps + integer(kind=i_def) :: ls_nstage, ls_substep + integer(kind=i_def) :: ls_stage, number_substeps + integer(kind=i_def) :: ls_number_substeps integer(kind=i_def) :: step integer(kind=i_def) :: splitting real(kind=r_def) :: dt_mol_substep @@ -102,7 +105,9 @@ module atl_mol_conservative_alg_mod type(field_type), allocatable :: rk_field(:) real(kind=r_def), allocatable :: rk_weights(:,:) type(transport_metadata_type), pointer :: transport_metadata + type(transport_metadata_type), pointer :: ls_transport_metadata type(transport_counter_type), pointer :: transport_counter + type(transport_counter_type), pointer :: ls_transport_counter type(transport_controller_type), pointer :: transport_controller type(transport_controller_type), pointer :: ls_transport_controller type(transport_controller_type), pointer :: pert_transport_controller @@ -122,11 +127,18 @@ module atl_mol_conservative_alg_mod transport_controller => tl_transport_controller%get_ls_wind_ls_rho_controller() ls_transport_controller => tl_transport_controller%get_ls_wind_pert_rho_controller() pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() + transport_metadata => pert_transport_controller%get_transport_metadata() + ls_transport_metadata => transport_controller%get_transport_metadata() + transport_counter => pert_transport_controller%get_transport_counter() + ls_transport_counter => transport_controller%get_transport_counter() + flux_precomputations => transport_controller%get_flux_precomputations() + ls_wind_flux_precomp => ls_transport_controller%get_flux_precomputations() pert_wind_flux_precomp => pert_transport_controller%get_flux_precomputations() + ls_wind_precomputations => ls_transport_controller%get_wind_precomputations() pert_wind_precomputations => pert_transport_controller%get_wind_precomputations() @@ -142,16 +154,27 @@ module atl_mol_conservative_alg_mod dt_mol_substep = ls_wind_precomputations%get_dt_mol_substep( & mesh%get_id(), direction, splitting, step & ) - div => get_directional_im3_div(mesh_id, direction) + div => get_directional_im3_div_r_tran(mesh_id, direction) call get_rk_transport_weights(nstage, rk_weights, runge_kutta_method) + + if (transport_efficiency) then + ls_nstage = 1 + ls_number_substeps = 1 + else + ls_nstage = nstage + ls_number_substeps = number_substeps + end if + allocate( rk_field(nstage) ) - allocate( stored_ls_field(number_substeps,nstage) ) + allocate( stored_ls_field(ls_number_substeps,ls_nstage) ) do stage = 1, nstage call rk_field(stage)%initialise(field%get_function_space()) - do substep = 1, number_substeps - call stored_ls_field(substep, stage)%initialise(field%get_function_space()) + end do + do ls_stage = 1, ls_nstage + do ls_substep = 1, ls_number_substeps + call stored_ls_field(ls_substep, ls_stage)%initialise(field%get_function_space()) end do end do call rhs_field%initialise( field%get_function_space() ) @@ -186,20 +209,24 @@ module atl_mol_conservative_alg_mod ! This calculates the values of ls_field. These are stored as an ! array for each substep and each RK stage to use in the perturbation ! ======================================================================== ! + + if ( transport_efficiency ) then + call invoke( setval_X( stored_ls_field(1,1), ls_field_np1 ) ) + else ! ------------------------------------------------------------------------ ! ! Start of substepping loop ! ------------------------------------------------------------------------ ! - ls_substep_loop: do substep = 1, number_substeps + ls_substep_loop: do substep = 1, ls_number_substeps - final_substep = ( substep == number_substeps ) + final_substep = ( substep == ls_number_substeps ) ! Reset field_n ready for the this substep call ls_field_np1%copy_field_properties(ls_field_n) call invoke( setval_X(ls_field_n, ls_field_np1) ) - ls_stage_loop: do stage = 1, nstage + ls_stage_loop: do stage = 1, ls_nstage - final_rk_stage = ( stage == nstage ) + final_rk_stage = ( stage == ls_nstage ) ! Store values for use in the perturbation code. This is done before ! advecting so as to pick up the correct linearisation state @@ -240,7 +267,7 @@ module atl_mol_conservative_alg_mod call advective_and_flux_alg(flux_step, adv_inc, rhs_field, ls_field_n, & ls_advecting_wind, direction, & - transport_metadata, final_rk_stage, & + ls_transport_metadata, final_rk_stage, & dt_mol_substep, do_flux, do_advective) if ( do_flux ) then @@ -283,9 +310,11 @@ module atl_mol_conservative_alg_mod ! increment, and may be adjusted to enforce min value or in blending zone call end_of_conservative_step_alg( & ls_field_np1, ls_field, ls_sum_flux, flux_precomputations, & - transport_counter, transport_metadata & + ls_transport_counter, ls_transport_metadata & ) + end if + ! ======================================================================== ! ! Perturbation ! @@ -323,6 +352,12 @@ module atl_mol_conservative_alg_mod ! ------------------------------------------------------------------------ ! substep_loop: do substep = number_substeps, 1, -1 + if (transport_efficiency) then + ls_substep = 1 + else + ls_substep = substep + end if + final_substep = ( substep == number_substeps ) ! ---------------------------------------------------------------------- ! @@ -339,6 +374,12 @@ module atl_mol_conservative_alg_mod stage_loop: do stage = nstage, 1, -1 + if (transport_efficiency) then + ls_stage = 1 + else + ls_stage = stage + end if + final_rk_stage = ( stage == nstage ) ! -------------------------------------------------------------------- ! @@ -387,7 +428,7 @@ module atl_mol_conservative_alg_mod call atl_advective_and_flux_alg( & flux_step_ls_wind, flux_step_pert_wind, adv_inc, & rhs_field, advecting_wind, & - stored_ls_field(substep,stage), ls_advecting_wind, & + stored_ls_field(ls_substep,ls_stage), ls_advecting_wind, & direction, transport_metadata, final_rk_stage, & dt_mol_substep, do_flux, do_advective & ) diff --git a/science/gungho/source/algorithm/transport/common/end_of_transport_step_alg_mod.x90 b/science/gungho/source/algorithm/transport/common/end_of_transport_step_alg_mod.x90 index 34ee74e25..ec2e5e80e 100644 --- a/science/gungho/source/algorithm/transport/common/end_of_transport_step_alg_mod.x90 +++ b/science/gungho/source/algorithm/transport/common/end_of_transport_step_alg_mod.x90 @@ -837,6 +837,7 @@ contains ! General case for any splitting else + call sum_hori_flux%initialise(w2h_fs) call sum_vert_flux%initialise(w2v_fs) diff --git a/science/gungho/source/algorithm/transport/common/transport_counter_mod.x90 b/science/gungho/source/algorithm/transport/common/transport_counter_mod.x90 index 6e7dbe6a6..b40c2cc84 100644 --- a/science/gungho/source/algorithm/transport/common/transport_counter_mod.x90 +++ b/science/gungho/source/algorithm/transport/common/transport_counter_mod.x90 @@ -40,6 +40,7 @@ module transport_counter_mod contains procedure, public :: initialise + procedure, public :: adj_initialise procedure, public :: get_field_n procedure, public :: set_field_n procedure, public :: get_flux @@ -51,6 +52,7 @@ module transport_counter_mod procedure, public :: get_num_split_steps_per_substep procedure, public :: get_dt_substep procedure, public :: get_outer_iteration + procedure, public :: dec_split_step_counter procedure, public :: inc_split_step_counter procedure, public :: inc_substep_counter procedure, public :: apply_cheap_update @@ -98,18 +100,70 @@ contains self%split_step_of_substep_counter = 1 self%substep_counter = 1 - if (allocated(self%flux)) deallocate(self%flux) - equation_form = transport_metadata%get_equation_form() ! Allocate arrays for dry fields if (equation_form == equation_form_conservative & - .or. equation_form == equation_form_consistent) then - allocate(self%flux(self%num_split_steps_per_substep-1)) + .or. equation_form == equation_form_consistent) then + if ( .not. allocated(self%flux)) then + allocate(self%flux(self%num_split_steps_per_substep-1)) + end if end if end subroutine initialise + !> @brief Initialises the transport_counter object, allocating flux fields + !! if necessary, ready for the adjoint model simulation + !> @details This is the same as initialise, but it prepares the code to be + !! run in reverse order, with the counter decreasing instead of + !! increasing. + !> @param[in] transport_metadata Object containing metadata describing + !! the options for transporting a variable + !> @param[in] dt_substep The time interval for a single substep + !> @param[in] num_substeps The number of total transport substeps to be + !! performed + !> @param[in] outer Optional, current iteration of the outer + !! loop of the semi-implicit time step + subroutine adj_initialise(self, transport_metadata, dt_substep, num_substeps, outer) + + implicit none + + class(transport_counter_type), intent(inout) :: self + type(transport_metadata_type), intent(in) :: transport_metadata + real(kind=r_tran), intent(in) :: dt_substep + integer(kind=i_def), intent(in) :: num_substeps + integer(kind=i_def), optional, intent(in) :: outer + + integer(kind=i_def) :: equation_form + + self%dt_substep = dt_substep + self%num_substeps = num_substeps + self%num_split_steps_per_substep = get_num_split_steps(transport_metadata%get_splitting()) + self%num_split_steps_per_whole_step = self%num_substeps * self%num_split_steps_per_substep + + if (present(outer)) then + self%outer = outer + else + self%outer = 0 + end if + + ! Initialise step counters to the last values (NOT 1) + self%split_step_of_whole_step_counter = self%num_split_steps_per_whole_step + self%split_step_of_substep_counter = self%num_split_steps_per_substep + self%substep_counter = num_substeps + + equation_form = transport_metadata%get_equation_form() + + ! Allocate arrays for dry fields + if (equation_form == equation_form_conservative & + .or. equation_form == equation_form_consistent) then + if ( .not. allocated(self%flux)) then + allocate(self%flux(self%num_split_steps_per_substep-1)) + end if + end if + + end subroutine adj_initialise + !> Public finalise method for the transport_counter subroutine finalise(self) @@ -314,6 +368,18 @@ contains ! ============================================================================ ! ! UTILITIES ! ============================================================================ ! + !> @brief Decreases the counter of transport split steps (for adjoint model) + subroutine dec_split_step_counter(self) + + implicit none + + class(transport_counter_type), target, intent(inout) :: self + + self%split_step_of_substep_counter = self%split_step_of_substep_counter - 1 + self%split_step_of_whole_step_counter = self%split_step_of_whole_step_counter - 1 + + end subroutine dec_split_step_counter + !> @brief Increments the counter of transport split steps subroutine inc_split_step_counter(self) diff --git a/science/linear/integration-test/runge_kutta/resources/runge_kutta_configuration.nml b/science/linear/integration-test/runge_kutta/resources/runge_kutta_configuration.nml index 9298f96f7..0296f5165 100644 --- a/science/linear/integration-test/runge_kutta/resources/runge_kutta_configuration.nml +++ b/science/linear/integration-test/runge_kutta/resources/runge_kutta_configuration.nml @@ -139,6 +139,7 @@ write_minmax_tseries=.false., &linear fixed_ls=.false. pert_option='random', +transport_efficiency=.false. / &logging run_log_level='info', diff --git a/science/linear/integration-test/semi_implicit/resources/semi_implicit_configuration.nml b/science/linear/integration-test/semi_implicit/resources/semi_implicit_configuration.nml index 2b208a451..902c03a53 100644 --- a/science/linear/integration-test/semi_implicit/resources/semi_implicit_configuration.nml +++ b/science/linear/integration-test/semi_implicit/resources/semi_implicit_configuration.nml @@ -140,6 +140,7 @@ write_minmax_tseries=.false., &linear fixed_ls=.false. pert_option='random', +transport_efficiency=.false. / &logging run_log_level='info', diff --git a/science/linear/rose-meta/lfric-linear/HEAD/rose-meta.conf b/science/linear/rose-meta/lfric-linear/HEAD/rose-meta.conf index 98a373e47..2e8b47047 100644 --- a/science/linear/rose-meta/lfric-linear/HEAD/rose-meta.conf +++ b/science/linear/rose-meta/lfric-linear/HEAD/rose-meta.conf @@ -4,7 +4,7 @@ import=lfric-gungho/HEAD compulsory=true description=Provides extra information for the tangent linear model. help= -ns=namelist/Linear/Initial +ns=namelist/Linear sort-key=Section-A07 title=Linear @@ -13,6 +13,7 @@ compulsory=true description=Do not update the linearisation state during the timestep. help=The linearisation state is the same for every outer and inner loop. !kind=default +ns=namelist/Linear/Efficiency type=logical [namelist:linear=l_stabilise_bl] @@ -20,7 +21,7 @@ compulsory=true description=?????? help=If true then turn on boundary layer stabilisation in linear model !kind=default -ns=namelist/Job/Timestepping/semi-implicit +ns=namelist/Linear/BoundaryLayer type=logical [namelist:linear=ls_read_w2h] @@ -30,7 +31,7 @@ description=For the linearisation state, read winds on a W2h space =from a file. help=Horizontal winds are read in on a W2h space from the ls file. If =false, then cell-centres (W3), co-located winds are assumed. -sort-key=Panel-A07a +ns=namelist/Linear/Data type=logical [namelist:linear=max_bl_stabilisation] @@ -39,7 +40,7 @@ description=?????? fail-if=this < 0.0 help=Maximum boundary layer stabilisation scalar in linear model !kind=default -ns=namelist/Job/Timestepping/semi-implicit +ns=namelist/Linear/BoundaryLayer range=0.0: type=real @@ -49,7 +50,7 @@ description=?????? fail-if=this < 1 ; help=Number of boundary layers levels to stabilise in linear model !kind=default -ns=namelist/Job/Timestepping/semi-implicit +ns=namelist/Linear/BoundaryLayer range=1: type=integer @@ -61,20 +62,26 @@ help=This can be an analytical perturbation, random data or from a file. =The analytical perturbation is based on the DCMIP gravity wave test. =The random data is generated with the Fortran intrinsic function 'random_number'. =The file contains data produced from LFRic with horizontal winds in W2. -ns=namelist/Linear/Initial -sort-key=Panel-A07a +ns=namelist/Linear/Data value-titles=analytic, =random, =file, =zero values='analytic', 'random', 'file', 'zero' +[namelist:linear=transport_efficiency] +compulsory=true +description=Do not update the linearisation state during the transport step. +help=The linearisation state remains fixed. +!kind=default +ns=namelist/Linear/Efficiency +type=logical + [namelist:validity_test] compulsory=true description=Provides settings for the tangent linear validity tests. help=The validity tests are performed using the integration test framework. ns=namelist/Linear/ValidityTest -sort-key=Section-A07 title=Validity Tests [namelist:validity_test=number_gamma_values] @@ -87,7 +94,6 @@ help=Gamma is the size of the perturbation used in the validity test. !kind=default ns=namelist/Linear/ValidityTest range=0: -sort-key=Panel-A07b type=integer [namelist:validity_test=update_ls_frequency] @@ -102,5 +108,4 @@ help=In the integration tests, when running over multiple timesteps, !kind=default ns=namelist/Linear/ValidityTest range=0: -sort-key=Panel-A07c type=integer diff --git a/science/linear/rose-meta/lfric-linear/versions.py b/science/linear/rose-meta/lfric-linear/versions.py index 152c043d0..90fa49f2e 100644 --- a/science/linear/rose-meta/lfric-linear/versions.py +++ b/science/linear/rose-meta/lfric-linear/versions.py @@ -20,14 +20,34 @@ def __repr__(self): """ Copy this template and complete to add your macro - class vnXX_txxx(MacroUpgrade): # Upgrade macro for by - BEFORE_TAG = "vnX.X" AFTER_TAG = "vnX.X_txxx" - def upgrade(self, config, meta_config=None): # Add settings return config, self.reports """ + + +class vn30_t108(MacroUpgrade): + """Upgrade macro for ticket #108 by Christine Johnson.""" + + BEFORE_TAG = "vn3.0" + AFTER_TAG = "vn3.0_t108" + + def upgrade(self, config, meta_config=None): + # Commands From: rose-meta/lfric-linear + fixed_ls = self.get_setting_value( + config, ["namelist:linear", "fixed_ls"] + ) + if ".true." in fixed_ls: + self.add_setting( + config, ["namelist:linear", "transport_efficiency"], ".true." + ) + else: + self.add_setting( + config, ["namelist:linear", "transport_efficiency"], ".false." + ) + + return config, self.reports diff --git a/science/linear/source/algorithm/transport/control/tl_split_transport_mod.x90 b/science/linear/source/algorithm/transport/control/tl_split_transport_mod.x90 new file mode 100644 index 000000000..779dd5b5a --- /dev/null +++ b/science/linear/source/algorithm/transport/control/tl_split_transport_mod.x90 @@ -0,0 +1,333 @@ +!----------------------------------------------------------------------------- +! (C) Crown copyright 2026 Met Office. All rights reserved. +! The file LICENCE, distributed with this code, contains details of the terms +! under which the code may be used. +!----------------------------------------------------------------------------- +!> @brief Tangent linear version of split_transport + +module tl_split_transport_mod + + use constants_mod, only: i_def, r_tran + use r_tran_field_mod, only: r_tran_field_type + use log_mod, only: log_event, LOG_LEVEL_ERROR + use split_transport_utils_mod, only: get_splitting_direction, & + get_splitting_fraction, & + get_num_split_steps + use split_transport_mod, only: split_transport_field + + ! Transport control infrastructure + use transport_controller_mod, only: transport_controller_type + use tl_transport_controller_mod, only: tl_transport_controller_type + use transport_counter_mod, only: transport_counter_type + use transport_enumerated_types_mod, only: direction_h, & + direction_v, & + split_method_null, & + split_method_mol, & + split_method_ffsl, & + split_method_sl, & + equation_form_conservative, & + equation_form_advective, & + equation_form_consistent + use transport_metadata_mod, only: transport_metadata_type + + ! Algorithms + use end_of_transport_step_alg_mod, only: end_of_conservative_step_alg, & + end_of_consistent_step_alg, & + end_of_advective_step_alg + use linear_config_mod, only: transport_efficiency + + implicit none + + private + + public :: tl_split_transport_control + public :: tl_split_transport_field + +contains + + !============================================================================= + !> @brief Controls vertical/horizontal split transport. + !> @details Manages the vertical/horizontal splitting of the split transport + !! scheme by choosing the splitting type and calling the individual + !! vertical and horizontal split components. + !> @param[in,out] field_np1 ACTIVE Field to return at end + !> @param[in] field_n ACTIVE Field at the start + !> @param[in] ls_field_n PASSIVE LS Field at the start + !> @param[in,out] tl_transport_controller + !! Encapsulating object containing the + !! transport counter and precomputations + subroutine tl_split_transport_control(field_np1, field_n, ls_field_n, & + tl_transport_controller) + + implicit none + + ! Arguments + type(r_tran_field_type), intent(inout) :: field_np1 + type(r_tran_field_type), target, intent(in) :: field_n + type(r_tran_field_type), target, intent(in) :: ls_field_n + type(tl_transport_controller_type), intent(inout) :: tl_transport_controller + + ! Internal variables + integer(kind=i_def) :: num_split_steps + integer(kind=i_def) :: split_step_count + type(r_tran_field_type), allocatable :: ls_field_np1(:) + type(r_tran_field_type), target :: field_tmp + type(r_tran_field_type), pointer :: field_ptr + type(transport_counter_type), pointer :: transport_counter + type(transport_counter_type), pointer :: ls_transport_counter + type(transport_metadata_type), pointer :: transport_metadata + type(transport_controller_type), pointer :: pert_transport_controller + type(transport_controller_type), pointer :: ls_transport_controller + real(kind=r_tran) :: dt_substep + integer(kind=i_def) :: num_substeps + + ! ------------------------------------------------------------------------ ! + ! NONLINEAR (LS) + ! ------------------------------------------------------------------------ ! + + ls_transport_controller => tl_transport_controller%get_ls_wind_ls_rho_controller() + transport_metadata => ls_transport_controller%get_transport_metadata() + ls_transport_counter => ls_transport_controller%get_transport_counter() + + ! Initialise fields + num_split_steps = get_num_split_steps(transport_metadata%get_splitting()) + allocate(ls_field_np1(num_split_steps + 1)) + do split_step_count = 1, num_split_steps + 1 + call field_n%copy_field_properties(ls_field_np1(split_step_count)) + end do + + if ( .not. transport_efficiency ) then + + ! When we have multiple split steps, we need an intermediate field for the + ! field at the start of each substep + ! field_ptr points to the field at the start of each split step + call invoke( setval_X(ls_field_np1(1), ls_field_n) ) + + if (num_split_steps > 1) then + call field_n%copy_field_properties(field_tmp) + call invoke( setval_X(field_tmp, ls_field_n) ) + field_ptr => field_tmp + else + field_ptr => ls_field_n + end if + + do split_step_count = 1, num_split_steps + call split_transport_field( ls_field_np1(split_step_count+1), & + field_ptr, & + ls_transport_controller & + ) + + if (split_step_count < num_split_steps) then + call invoke( setval_X(field_tmp, ls_field_np1(split_step_count+1) ) ) + + ! Increment split step counter + call ls_transport_counter%inc_split_step_counter() + end if + end do + + else + ! Assume that the linearisation state does not evolve during the split transport + do split_step_count = 1, num_split_steps + 1 + call invoke( setval_X(ls_field_np1(split_step_count), ls_field_n) ) + end do + end if + + dt_substep = ls_transport_counter%get_dt_substep() + num_substeps = ls_transport_counter%get_num_substeps() + call ls_transport_counter%initialise(transport_metadata, dt_substep, num_substeps ) + + ! ------------------------------------------------------------------------ ! + ! LINEAR + ! ------------------------------------------------------------------------ ! + + pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() + transport_metadata => pert_transport_controller%get_transport_metadata() + transport_counter => pert_transport_controller%get_transport_counter() + + ! Initialise fields + num_split_steps = get_num_split_steps(transport_metadata%get_splitting()) + + ! When we have multiple split steps, we need an intermediate field for the + ! field at the start of each substeps + ! field_ptr points to the field at the start of each split step + if (num_split_steps > 1) then + call field_n%copy_field_properties(field_tmp) + call invoke( setval_X(field_tmp, field_n) ) + field_ptr => field_tmp + else + field_ptr => field_n + end if + + do split_step_count = 1, num_split_steps + + call tl_split_transport_field( & + field_np1, field_ptr, ls_field_np1(split_step_count), & + tl_transport_controller & + ) + + if (split_step_count < num_split_steps) then + call invoke( setval_X(field_tmp, field_np1 ) ) + + ! Increment split step counter + call transport_counter%inc_split_step_counter() + call ls_transport_counter%inc_split_step_counter() + end if + + end do + + deallocate( ls_field_np1 ) + + end subroutine tl_split_transport_control + + !============================================================================= + !> @brief Does either vertical or horizontal transport of a field. + !> @details Performs a vertical or horizontal transport step, solving the + !! transport equation for a (multidata) field. + !> @param[in,out] field_np1 ACTIVE Field to return at end + !> @param[in] field_n ACTIVE Field at the start + !> @param[in] ls_field_n PASSIVE Field at the start + !> @param[in,out] tl_transport_controller + !! Encapsulating object containing the + !! transport counter and precomputations + subroutine tl_split_transport_field(field_np1, field_n, ls_field_n, tl_transport_controller) + + use tl_mol_conservative_alg_mod, only: tl_mol_conservative_alg + use tl_mol_advective_alg_mod, only: tl_mol_advective_alg + + implicit none + + ! Arguments + type(r_tran_field_type), intent(inout) :: field_np1 + type(r_tran_field_type), intent(in) :: field_n + type(r_tran_field_type), intent(in) :: ls_field_n + type(tl_transport_controller_type), intent(inout) :: tl_transport_controller + + ! Internal variables + integer(kind=i_def) :: method, direction + type(transport_counter_type), pointer :: transport_counter + type(transport_metadata_type), pointer :: transport_metadata + type(transport_controller_type), pointer :: pert_transport_controller + + pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() + transport_counter => pert_transport_controller%get_transport_counter() + transport_metadata => pert_transport_controller%get_transport_metadata() + + ! -------------------------------------------------------------------------! + ! Set up method based on direction + ! -------------------------------------------------------------------------! + direction = get_splitting_direction( & + transport_metadata%get_splitting(), & + transport_counter%get_split_step_of_substep_counter() & + ) + + select case ( direction ) + case ( direction_h ) + method = transport_metadata%get_horizontal_method() + case ( direction_v ) + method = transport_metadata%get_vertical_method() + case default + call log_event('Split transport direction not recognised', LOG_LEVEL_ERROR) + end select + + ! -------------------------------------------------------------------------! + ! Choose method, and then choose equation + ! -------------------------------------------------------------------------! + select case ( method ) + + ! -------------------------------------------------------------------------! + ! Null step + ! -------------------------------------------------------------------------! + case ( split_method_null ) + call log_event( & + 'TL: Split method null not coded', & + LOG_LEVEL_ERROR & + ) + + ! -------------------------------------------------------------------------! + ! Method of Lines step + ! -------------------------------------------------------------------------! + case ( split_method_mol ) + ! Choose form of transport equation + select case ( transport_metadata%get_equation_form() ) + case ( equation_form_conservative ) + call tl_mol_conservative_alg( & + field_np1, field_n, ls_field_n, tl_transport_controller & + ) + case ( equation_form_advective ) + call tl_mol_advective_alg( & + field_np1, field_n, ls_field_n, tl_transport_controller & + ) + case ( equation_form_consistent ) + call log_event( & + 'TL: MoL Consistent not coded yet', & + LOG_LEVEL_ERROR & + ) + case default + call log_event( & + 'Trying to solve unrecognised form of transport equation', & + LOG_LEVEL_ERROR & + ) + + end select + + ! -------------------------------------------------------------------------! + ! Flux-Form Semi-Lagrangian step + ! -------------------------------------------------------------------------! + case ( split_method_ffsl ) + ! All equation forms have the same control method + call log_event( & + 'TL: FFSL not coded yet', & + LOG_LEVEL_ERROR & + ) + ! -------------------------------------------------------------------------! + ! Semi-Lagrangian step + ! -------------------------------------------------------------------------! + case ( split_method_sl ) + ! Choose direction + select case ( direction ) + + case ( direction_h ) + ! Horizontal SL only for advective form + if ( transport_metadata%get_equation_form() /= equation_form_advective ) then + call log_event( & + 'Horizontal semi-Lagrangian is only for advective form', & + LOG_LEVEL_ERROR & + ) + end if + call log_event( & + 'TL: Horizontal SL not coded yet', & + LOG_LEVEL_ERROR & + ) + + case ( direction_v ) + ! Choose form of transport equation for vertical + select case ( transport_metadata%get_equation_form() ) + + case ( equation_form_conservative ) + call log_event( & + 'TL: Vertical SL conservative not coded yet', & + LOG_LEVEL_ERROR & + ) + case ( equation_form_advective ) + call log_event( & + 'TL: Vertical SL advective not coded yet', & + LOG_LEVEL_ERROR & + ) + case default + call log_event( & + 'Trying to solve unrecognised form of transport equation', & + LOG_LEVEL_ERROR & + ) + end select + end select + + case default + call log_event( & + 'Trying to transport with unrecognised scheme', & + LOG_LEVEL_ERROR & + ) + end select + + end subroutine tl_split_transport_field + +end module tl_split_transport_mod diff --git a/science/linear/source/algorithm/transport/control/tl_transport_field_mod.f90 b/science/linear/source/algorithm/transport/control/tl_transport_field_mod.f90 index 17016206a..5fd0a5b02 100644 --- a/science/linear/source/algorithm/transport/control/tl_transport_field_mod.f90 +++ b/science/linear/source/algorithm/transport/control/tl_transport_field_mod.f90 @@ -21,6 +21,7 @@ module tl_transport_field_mod equation_form_advective use tl_mol_conservative_alg_mod, only: tl_mol_conservative_alg use tl_mol_advective_alg_mod, only: tl_mol_advective_alg + use tl_split_transport_mod, only: tl_split_transport_control use transport_controller_mod, only: transport_controller_type use tl_transport_controller_mod, only: tl_transport_controller_type use transport_counter_mod, only: transport_counter_type @@ -121,10 +122,8 @@ subroutine tl_transport_field(field_np1, field_n, ls_field_n, & ! Some split horizontal/vertical transport scheme ! -------------------------------------------------------------------------! case ( scheme_split ) - call log_event( & - 'Split transport not implemented for tangent-linear app', & - LOG_LEVEL_ERROR & - ) + call tl_split_transport_control(field_np1, field_n, ls_field_n, & + tl_transport_controller) case default call log_event( & diff --git a/science/linear/source/algorithm/transport/mol/tl_mol_advective_alg_mod.x90 b/science/linear/source/algorithm/transport/mol/tl_mol_advective_alg_mod.x90 index a1a34919f..e85d71036 100644 --- a/science/linear/source/algorithm/transport/mol/tl_mol_advective_alg_mod.x90 +++ b/science/linear/source/algorithm/transport/mol/tl_mol_advective_alg_mod.x90 @@ -28,6 +28,7 @@ module tl_mol_advective_alg_mod ! Configuration use transport_config_mod, only: runge_kutta_method + use linear_config_mod, only: transport_efficiency implicit none @@ -58,7 +59,9 @@ contains ! Internal variables integer(kind=i_def) :: stage, s, direction integer(kind=i_def) :: nstage, substep + integer(kind=i_def) :: ls_nstage, ls_substep integer(kind=i_def) :: number_substeps + integer(kind=i_def) :: ls_stage, ls_number_substeps integer(kind=i_def) :: step integer(kind=i_def) :: splitting real(kind=r_def) :: dt_mol_substep @@ -73,7 +76,10 @@ contains type(field_type), allocatable :: rk_field(:) real(kind=r_def), allocatable :: rk_weights(:,:) type(transport_metadata_type), pointer :: transport_metadata + type(transport_metadata_type), pointer :: ls_transport_metadata type(transport_counter_type), pointer :: transport_counter + type(transport_counter_type), pointer :: ls_transport_counter + type(transport_controller_type), pointer :: transport_controller type(transport_controller_type), pointer :: ls_transport_controller type(transport_controller_type), pointer :: pert_transport_controller type(wind_precomputations_type), pointer :: ls_wind_precomputations @@ -83,10 +89,17 @@ contains ! Extract transport objects and initialise temporary fields ! ------------------------------------------------------------------------ ! mesh => field%get_mesh() + + transport_controller => tl_transport_controller%get_ls_wind_ls_rho_controller() ls_transport_controller => tl_transport_controller%get_ls_wind_pert_rho_controller() pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() + transport_metadata => pert_transport_controller%get_transport_metadata() + ls_transport_metadata => transport_controller%get_transport_metadata() + transport_counter => pert_transport_controller%get_transport_counter() + ls_transport_counter => transport_controller%get_transport_counter() + ls_wind_precomputations => ls_transport_controller%get_wind_precomputations() pert_wind_precomputations => pert_transport_controller%get_wind_precomputations() @@ -107,13 +120,23 @@ contains call field_n%initialise(field%get_function_space()) call get_rk_transport_weights(nstage, rk_weights, runge_kutta_method) + + if (transport_efficiency) then + ls_number_substeps = 1 + ls_nstage =1 + else + ls_number_substeps = number_substeps + ls_nstage = nstage + end if allocate( rk_field(nstage) ) - allocate( stored_ls_field(number_substeps, nstage) ) + allocate( stored_ls_field(ls_number_substeps, ls_nstage) ) do stage = 1, nstage call rk_field(stage)%initialise(field%get_function_space()) - do substep = 1, number_substeps - call stored_ls_field(substep, stage)%initialise(field%get_function_space()) + end do + do ls_stage = 1, ls_nstage + do ls_substep = 1, ls_number_substeps + call stored_ls_field(ls_substep, ls_stage)%initialise(field%get_function_space()) end do end do call rhs_field%initialise( field%get_function_space() ) @@ -130,19 +153,23 @@ contains ! array for each substep and each RK stage to use in the perturbation !--------------------------------------------------------------------------! + if (transport_efficiency) then + call invoke( setval_X( stored_ls_field(1,1), ls_field_np1)) + else + ! Perform the number of rk-stages and substeps required - do substep = 1, number_substeps + do substep = 1, ls_number_substeps ! Reset field_n ready for the this substep call ls_field_np1%copy_field_properties(ls_field_n) call invoke( setval_X(ls_field_n, ls_field_np1) ) - do stage = 1, nstage + do stage = 1, ls_nstage ! Store values for use in the perturbation code. This is done before ! advecting so as to pick up the correct linearisation state call invoke( setval_X(stored_ls_field(substep, stage), ls_field_np1) ) - final_rk_stage = ( stage == nstage ) + final_rk_stage = ( stage == ls_nstage ) ! Compute the field for this stage: ! rhs_field = sum(s=1,stage): a(stage,s)*field^(s) @@ -153,9 +180,9 @@ contains end do ! Compute update: rhs = u.grad(rhs_field) - call advective_and_flux_alg(dummy, rhs, rhs_field, ls_field_n, & - ls_advecting_wind, direction, & - transport_metadata, final_rk_stage, & + call advective_and_flux_alg(dummy, rhs, rhs_field, ls_field_n, & + ls_advecting_wind, direction, & + ls_transport_metadata, final_rk_stage, & dt_mol_substep, .false., .true. ) ! Update field: f = f^n - dt*rhs @@ -165,10 +192,12 @@ contains end do ! End of step: if necessary enforce min val and overwrite in blending zone - call end_of_advective_step_alg( & - ls_field_np1, ls_field_n, transport_counter, transport_metadata & + call end_of_advective_step_alg( & + ls_field_np1, ls_field_n, ls_transport_counter, ls_transport_metadata & ) + end if + ! -------------------------------------------------------------------------- ! Perturbation ! @@ -177,15 +206,28 @@ contains ! linearisation stage stage !--------------------------------------------------------------------------- + ! Perform the number of rk-stages and substeps required do substep = 1, number_substeps + if (transport_efficiency) then + ls_substep = 1 + else + ls_substep = substep + end if + ! Reset field_n ready for the this substep call field_np1%copy_field_properties(field_n) call invoke( setval_X(field_n, field_np1) ) do stage = 1, nstage + if (transport_efficiency) then + ls_stage = 1 + else + ls_stage = stage + end if + final_rk_stage = ( stage == nstage ) ! Compute the field for this stage: @@ -199,7 +241,7 @@ contains ! Compute update: rhs = u.grad(rhs_field) call tl_advective_and_flux_alg( & dummy, dummy, rhs, rhs_field, advecting_wind, & - stored_ls_field(substep, stage), ls_advecting_wind, & + stored_ls_field(ls_substep, ls_stage), ls_advecting_wind, & direction, transport_metadata, final_rk_stage, dt_mol_substep, & .false., .true. & ) diff --git a/science/linear/source/algorithm/transport/mol/tl_mol_conservative_alg_mod.x90 b/science/linear/source/algorithm/transport/mol/tl_mol_conservative_alg_mod.x90 index b4563234e..e9ec55497 100644 --- a/science/linear/source/algorithm/transport/mol/tl_mol_conservative_alg_mod.x90 +++ b/science/linear/source/algorithm/transport/mol/tl_mol_conservative_alg_mod.x90 @@ -27,7 +27,7 @@ module tl_mol_conservative_alg_mod use runge_kutta_init_mod, only: get_rk_transport_weights use split_transport_utils_mod, only: get_num_split_steps, & get_splitting_direction - use transport_constants_mod, only: get_directional_im3_div + use transport_constants_mod, only: get_directional_im3_div_r_tran use transport_controller_mod, only: transport_controller_type use tl_transport_controller_mod, only: tl_transport_controller_type use transport_counter_mod, only: transport_counter_type @@ -43,6 +43,7 @@ module tl_mol_conservative_alg_mod use boundaries_config_mod, only: limited_area use base_mesh_config_mod, only: topology, topology_non_periodic use io_config_mod, only: subroutine_timers + use linear_config_mod, only: transport_efficiency use transport_config_mod, only: runge_kutta_method, & dry_field_name, & operators, & @@ -79,6 +80,8 @@ module tl_mol_conservative_alg_mod integer(kind=i_def) :: stage, s, direction integer(kind=i_def) :: nstage, substep integer(kind=i_def) :: number_substeps + integer(kind=i_def) :: ls_nstage, ls_substep + integer(kind=i_def) :: ls_stage, ls_number_substeps integer(kind=i_def) :: step integer(kind=i_def) :: splitting real(kind=r_def) :: dt_mol_substep @@ -106,7 +109,9 @@ module tl_mol_conservative_alg_mod type(field_type), allocatable :: rk_field(:) real(kind=r_def), allocatable :: rk_weights(:,:) type(transport_metadata_type), pointer :: transport_metadata + type(transport_metadata_type), pointer :: ls_transport_metadata type(transport_counter_type), pointer :: transport_counter + type(transport_counter_type), pointer :: ls_transport_counter type(transport_controller_type), pointer :: transport_controller type(transport_controller_type), pointer :: ls_transport_controller type(transport_controller_type), pointer :: pert_transport_controller @@ -126,11 +131,18 @@ module tl_mol_conservative_alg_mod transport_controller => tl_transport_controller%get_ls_wind_ls_rho_controller() ls_transport_controller => tl_transport_controller%get_ls_wind_pert_rho_controller() pert_transport_controller => tl_transport_controller%get_pert_wind_ls_rho_controller() + transport_metadata => pert_transport_controller%get_transport_metadata() + ls_transport_metadata => transport_controller%get_transport_metadata() + transport_counter => pert_transport_controller%get_transport_counter() + ls_transport_counter => transport_controller%get_transport_counter() + flux_precomputations => transport_controller%get_flux_precomputations() + ls_wind_flux_precomp => ls_transport_controller%get_flux_precomputations() pert_wind_flux_precomp => pert_transport_controller%get_flux_precomputations() + ls_wind_precomputations => ls_transport_controller%get_wind_precomputations() pert_wind_precomputations => pert_transport_controller%get_wind_precomputations() @@ -146,16 +158,27 @@ module tl_mol_conservative_alg_mod dt_mol_substep = ls_wind_precomputations%get_dt_mol_substep( & mesh%get_id(), direction, splitting, step & ) - div => get_directional_im3_div(mesh_id, direction) + div => get_directional_im3_div_r_tran(mesh_id, direction) call get_rk_transport_weights(nstage, rk_weights, runge_kutta_method) + + if (transport_efficiency) then + ls_number_substeps = 1 + ls_nstage = 1 + else + ls_number_substeps = number_substeps + ls_nstage = nstage + end if + allocate( rk_field(nstage) ) - allocate( stored_ls_field(number_substeps,nstage) ) + allocate( stored_ls_field(ls_number_substeps,ls_nstage) ) do stage = 1, nstage call rk_field(stage)%initialise(field%get_function_space()) - do substep = 1, number_substeps - call stored_ls_field(substep, stage)%initialise(field%get_function_space()) + end do + do ls_stage = 1, ls_nstage + do ls_substep = 1, ls_number_substeps + call stored_ls_field(ls_substep, ls_stage)%initialise(field%get_function_space()) end do end do call rhs_field%initialise( field%get_function_space() ) @@ -194,20 +217,23 @@ module tl_mol_conservative_alg_mod ! array for each substep and each RK stage to use in the perturbation ! ======================================================================== ! + if (transport_efficiency) then + call invoke( setval_X( stored_ls_field(1,1), ls_field_np1 ) ) + else ! ------------------------------------------------------------------------ ! ! Start of substepping loop ! ------------------------------------------------------------------------ ! - ls_substep_loop: do substep = 1, number_substeps + ls_substep_loop: do substep = 1, ls_number_substeps - final_substep = ( substep == number_substeps ) + final_substep = ( substep == ls_number_substeps ) ! Reset field_n ready for the this substep call ls_field_np1%copy_field_properties(ls_field_n) call invoke( setval_X(ls_field_n, ls_field_np1) ) - ls_stage_loop: do stage = 1, nstage + ls_stage_loop: do stage = 1, ls_nstage - final_rk_stage = ( stage == nstage ) + final_rk_stage = ( stage == ls_nstage ) ! Store values for use in the perturbation code. This is done before ! advecting so as to pick up the correct linearisation state @@ -248,7 +274,7 @@ module tl_mol_conservative_alg_mod call advective_and_flux_alg(flux_step, adv_inc, rhs_field, ls_field_n, & ls_advecting_wind, direction, & - transport_metadata, final_rk_stage, & + ls_transport_metadata, final_rk_stage, & dt_mol_substep, do_flux, do_advective) if ( do_flux ) then @@ -291,9 +317,11 @@ module tl_mol_conservative_alg_mod ! increment, and may be adjusted to enforce min value or in blending zone call end_of_conservative_step_alg( & ls_field_np1, ls_field, ls_sum_flux, flux_precomputations, & - transport_counter, transport_metadata & + ls_transport_counter, ls_transport_metadata & ) + end if + ! ======================================================================== ! ! Perturbation ! @@ -307,6 +335,12 @@ module tl_mol_conservative_alg_mod ! ------------------------------------------------------------------------ ! substep_loop: do substep = 1, number_substeps + if (transport_efficiency) then + ls_substep = 1 + else + ls_substep = substep + end if + final_substep = ( substep == number_substeps ) ! Reset field_n ready for the this substep @@ -315,6 +349,12 @@ module tl_mol_conservative_alg_mod stage_loop: do stage = 1, nstage + if (transport_efficiency) then + ls_stage = 1 + else + ls_stage = stage + end if + final_rk_stage = ( stage == nstage ) ! Compute the field for this stage: @@ -352,7 +392,7 @@ module tl_mol_conservative_alg_mod call tl_advective_and_flux_alg( & flux_step_ls_wind, flux_step_pert_wind, adv_inc, & rhs_field, advecting_wind, & - stored_ls_field(substep,stage), ls_advecting_wind, & + stored_ls_field(ls_substep,ls_stage), ls_advecting_wind, & direction, transport_metadata, final_rk_stage, & dt_mol_substep, do_flux, do_advective & )