From ce546c7544c281e3b597a3f71f087ffd2ab77721 Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Wed, 24 Dec 2025 07:45:15 -0700 Subject: [PATCH] pass diag fields (HS, T01, T0M1, THM) to cpl --- model/src/wav_comp_nuopc.F90 | 24 +++++++++++ model/src/wav_import_export.F90 | 75 ++++++++++++++++++++++++++++++++- 2 files changed, 98 insertions(+), 1 deletion(-) diff --git a/model/src/wav_comp_nuopc.F90 b/model/src/wav_comp_nuopc.F90 index bb42fde8f0..c08b179db0 100644 --- a/model/src/wav_comp_nuopc.F90 +++ b/model/src/wav_comp_nuopc.F90 @@ -920,6 +920,10 @@ subroutine DataInitialize(gcomp, rc) real(r8), pointer :: sw_lasl(:) real(r8), pointer :: sw_ustokes(:) real(r8), pointer :: sw_vstokes(:) + real(r8), pointer :: Sw_Hs(:) + real(r8), pointer :: Sw_t01(:) + real(r8), pointer :: Sw_t0m1(:) + real(r8), pointer :: Sw_thm(:) real(r8), pointer :: wave_elevation_spectrum(:,:) character(len=*),parameter :: subname = '(wav_comp_nuopc:DataInitialize)' ! ------------------------------------------------------------------- @@ -967,6 +971,26 @@ subroutine DataInitialize(gcomp, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return wave_elevation_spectrum(:,:) = 0. endif + if (state_fldchk(exportState, 'Sw_Hs')) then + call state_getfldptr(exportState, 'Sw_Hs', Sw_Hs, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + Sw_Hs (:) = 0. + endif + if (state_fldchk(exportState, 'Sw_t01')) then + call state_getfldptr(exportState, 'Sw_t01', Sw_t01, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + Sw_t01 (:) = 0. + endif + if (state_fldchk(exportState, 'Sw_t0m1')) then + call state_getfldptr(exportState, 'Sw_t0m1', Sw_t0m1, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + Sw_t0m1 (:) = 0. + endif + if (state_fldchk(exportState, 'Sw_thm')) then + call state_getfldptr(exportState, 'Sw_thm', Sw_thm, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + Sw_thm (:) = 0. + endif if (.not. unstr_mesh) then ! Set global grid size scalars in export state diff --git a/model/src/wav_import_export.F90 b/model/src/wav_import_export.F90 index 4b09c74457..5e0a106ee6 100644 --- a/model/src/wav_import_export.F90 +++ b/model/src/wav_import_export.F90 @@ -148,6 +148,10 @@ subroutine advertise_fields(importState, ExportState, flds_scalar_name, rc) end if call fldlist_add(fldsFrWav_num, fldsFrWav, 'Sw_pstokes_x', ungridded_lbound=1, ungridded_ubound=3) call fldlist_add(fldsFrWav_num, fldsFrWav, 'Sw_pstokes_y', ungridded_lbound=1, ungridded_ubound=3) + call fldlist_add(fldsFrWav_num, fldsFrWav, 'Sw_Hs') + call fldlist_add(fldsFrWav_num, fldsFrWav, 'Sw_t01') + call fldlist_add(fldsFrWav_num, fldsFrWav, 'Sw_t0m1') + call fldlist_add(fldsFrWav_num, fldsFrWav, 'Sw_thm') ! AA TODO: In the above fldlist_add calls, we are passing hardcoded ungridded_ubound values (3) because, USSPF(2) ! is not initialized yet. It is set during w3init which gets called at a later phase (realize). A permanent solution @@ -589,7 +593,7 @@ subroutine export_fields (gcomp, rc) use w3iogomd , only : CALC_U3STOKES #ifdef W3_CESMCOUPLED use w3wdatmd , only : ASF, UST - use w3adatmd , only : USSHX, USSHY, UD, HS + use w3adatmd , only : USSHX, USSHY, UD, HS, T01, T0M1, THM use w3idatmd , only : HSL #else use wmmdatmd , only : mdse, mdst, wmsetm @@ -622,6 +626,10 @@ subroutine export_fields (gcomp, rc) !real(r8), pointer :: sw_lasl(:) real(r8), pointer :: sw_ustokes(:) real(r8), pointer :: sw_vstokes(:) + real(r8), pointer :: sw_hs(:) + real(r8), pointer :: sw_t01(:) + real(r8), pointer :: sw_t0m1(:) + real(r8), pointer :: sw_thm(:) ! d2 is location, d1 is frequency - nwav_elev_spectrum frequencies will be used real(r8), pointer :: wave_elevation_spectrum(:,:) @@ -741,6 +749,71 @@ subroutine export_fields (gcomp, rc) call CalcRoughl(z0rlen) endif + if (state_fldchk(exportState, 'Sw_Hs')) then + call state_getfldptr(exportState, 'Sw_Hs', Sw_Hs, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + Sw_Hs(:) = fillvalue + do jsea=1, nseal_cpl + call init_get_isea(isea, jsea) + ix = mapsf(isea,1) + iy = mapsf(isea,2) + if (mapsta(iy,ix) == 1) then + Sw_Hs(jsea) = HS(jsea) + else + Sw_Hs(jsea) = 0. + endif + enddo + end if + + if (state_fldchk(exportState, 'Sw_t01')) then + call state_getfldptr(exportState, 'Sw_t01', sw_t01, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + sw_t01(:) = fillvalue + do jsea=1, nseal_cpl + call init_get_isea(isea, jsea) + ix = mapsf(isea,1) + iy = mapsf(isea,2) + if (mapsta(iy,ix) == 1) then + sw_t01(jsea) = t01(jsea) + else + sw_t01(jsea) = 0. + endif + enddo + end if + + if (state_fldchk(exportState, 'Sw_t0m1')) then + call state_getfldptr(exportState, 'Sw_t0m1', sw_t0m1, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + sw_t0m1(:) = fillvalue + do jsea=1, nseal_cpl + call init_get_isea(isea, jsea) + ix = mapsf(isea,1) + iy = mapsf(isea,2) + if (mapsta(iy,ix) == 1) then + sw_t0m1(jsea) = t0m1(jsea) + else + sw_t0m1(jsea) = 0. + endif + enddo + end if + + if (state_fldchk(exportState, 'Sw_thm')) then + call state_getfldptr(exportState, 'Sw_thm', sw_thm, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + sw_thm(:) = fillvalue + do jsea=1, nseal_cpl + call init_get_isea(isea, jsea) + ix = mapsf(isea,1) + iy = mapsf(isea,2) + if (mapsta(iy,ix) == 1) then + sw_thm(jsea) = thm(jsea) + else + sw_thm(jsea) = 0. + endif + enddo + end if + + if ( state_fldchk(exportState, 'wbcuru') .and. & state_fldchk(exportState, 'wbcurv') .and. & state_fldchk(exportState, 'wbcurp')) then