From 6ed235904f9e2c61d1caafc84f4d3386a65c930a Mon Sep 17 00:00:00 2001 From: JFLemieux73 Date: Fri, 2 Mar 2018 19:16:27 +0000 Subject: [PATCH] Added pressure, modified norm of principal stresses and made small modifs to basal stress following Till's comments --- cicecore/cicedynB/analysis/ice_history.F90 | 22 ++++++++--- .../cicedynB/analysis/ice_history_shared.F90 | 5 ++- cicecore/cicedynB/dynamics/ice_dyn_eap.F90 | 12 ++---- cicecore/cicedynB/dynamics/ice_dyn_evp.F90 | 11 ++---- cicecore/cicedynB/dynamics/ice_dyn_shared.F90 | 38 +++++++++++-------- cicecore/cicedynB/general/ice_flux.F90 | 7 ++-- .../io/io_binary/ice_history_write.F90 | 2 +- .../io/io_netcdf/ice_history_write.F90 | 5 ++- .../io/io_pio/ice_history_write.F90 | 5 ++- configuration/scripts/ice_in | 1 + 10 files changed, 59 insertions(+), 49 deletions(-) diff --git a/cicecore/cicedynB/analysis/ice_history.F90 b/cicecore/cicedynB/analysis/ice_history.F90 index 32dbe5823..85567ff7e 100644 --- a/cicecore/cicedynB/analysis/ice_history.F90 +++ b/cicecore/cicedynB/analysis/ice_history.F90 @@ -5,7 +5,7 @@ ! ! The following variables are currently hard-wired as snapshots ! (instantaneous rather than time-averages): -! divu, shear, sig1, sig2, trsig, mlt_onset, frz_onset, hisnap, aisnap +! divu, shear, sig1, sig2, sigI, trsig, mlt_onset, frz_onset, hisnap, aisnap ! ! Options for histfreq: '1','h','d','m','y','x', where x means that ! output stream will not be used (recommended for efficiency). @@ -300,6 +300,7 @@ subroutine init_hist (dt) call broadcast_scalar (f_shear, master_task) call broadcast_scalar (f_sig1, master_task) call broadcast_scalar (f_sig2, master_task) + call broadcast_scalar (f_sigI, master_task) call broadcast_scalar (f_dvidtt, master_task) call broadcast_scalar (f_dvidtd, master_task) call broadcast_scalar (f_daidtt, master_task) @@ -776,6 +777,11 @@ subroutine init_hist (dt) "norm. principal stress 2", & "sig2 is instantaneous", c1, c0, & ns1, f_sig2) + + call define_hist_field(n_sigI,"sigI","1",ustr2D, ucstr, & + "ice pressure", & + "sigI is instantaneous", c1, c0, & + ns1, f_sigI) call define_hist_field(n_dvidtt,"dvidtt","cm/day",tstr2D, tcstr, & "volume tendency thermo", & @@ -1199,11 +1205,11 @@ subroutine accum_hist (dt) fhocn, fhocn_ai, uatm, vatm, fbot, & fswthru_ai, strairx, strairy, strtltx, strtlty, strintx, strinty, & taubx, tauby, strocnx, strocny, fm, daidtt, dvidtt, daidtd, dvidtd, fsurf, & - fcondtop, fsurfn, fcondtopn, flatn, fsensn, albcnt, prs_sig, & + fcondtop, fsurfn, fcondtopn, flatn, fsensn, albcnt, & stressp_1, stressm_1, stress12_1, & stressp_2, stressm_2, stress12_2, & stressp_3, stressm_3, stress12_3, & - stressp_4, stressm_4, stress12_4, sig1, sig2, & + stressp_4, stressm_4, stress12_4, sig1, sig2, sigI, & mlt_onset, frz_onset, dagedtt, dagedtd, fswint_ai, keffn_top, & snowfrac, alvdr_ai, alvdf_ai, alidr_ai, alidf_ai use ice_arrays_column, only: snowfracn @@ -1922,9 +1928,10 @@ subroutine accum_hist (dt) stressp_1 (:,:,iblk), & stressm_1 (:,:,iblk), & stress12_1(:,:,iblk), & - prs_sig (:,:,iblk), & + strength (:,:,iblk), & sig1 (:,:,iblk), & - sig2 (:,:,iblk)) + sig2 (:,:,iblk), & + sigI (:,:,iblk)) do j = jlo, jhi do i = ilo, ihi @@ -1933,6 +1940,7 @@ subroutine accum_hist (dt) if (n_shear (ns) /= 0) a2D(i,j,n_shear(ns), iblk) = spval if (n_sig1 (ns) /= 0) a2D(i,j,n_sig1(ns), iblk) = spval if (n_sig2 (ns) /= 0) a2D(i,j,n_sig2(ns), iblk) = spval + if (n_sigI (ns) /= 0) a2D(i,j,n_sigI(ns), iblk) = spval if (n_mlt_onset(ns) /= 0) a2D(i,j,n_mlt_onset(ns),iblk) = spval if (n_frz_onset(ns) /= 0) a2D(i,j,n_frz_onset(ns),iblk) = spval if (n_hisnap (ns) /= 0) a2D(i,j,n_hisnap(ns), iblk) = spval @@ -1961,6 +1969,8 @@ subroutine accum_hist (dt) sig1 (i,j,iblk)*avail_hist_fields(n_sig1(ns))%cona if (n_sig2 (ns) /= 0) a2D(i,j,n_sig2(ns),iblk) = & sig2 (i,j,iblk)*avail_hist_fields(n_sig2(ns))%cona + if (n_sigI (ns) /= 0) a2D(i,j,n_sigI(ns),iblk) = & + sigI (i,j,iblk)*avail_hist_fields(n_sigI(ns))%cona if (n_mlt_onset(ns) /= 0) a2D(i,j,n_mlt_onset(ns),iblk) = & mlt_onset(i,j,iblk) if (n_frz_onset(ns) /= 0) a2D(i,j,n_frz_onset(ns),iblk) = & @@ -1972,7 +1982,7 @@ subroutine accum_hist (dt) if (kdyn == 2) then ! for EAP dynamics different time of output if (n_trsig (ns) /= 0) a2D(i,j,n_trsig(ns),iblk ) = & - prs_sig(i,j,iblk) + strength(i,j,iblk) else if (n_trsig (ns) /= 0) a2D(i,j,n_trsig(ns),iblk ) = & p25*(stressp_1(i,j,iblk) & diff --git a/cicecore/cicedynB/analysis/ice_history_shared.F90 b/cicecore/cicedynB/analysis/ice_history_shared.F90 index f98f10044..633fc314e 100644 --- a/cicecore/cicedynB/analysis/ice_history_shared.F90 +++ b/cicecore/cicedynB/analysis/ice_history_shared.F90 @@ -5,7 +5,7 @@ ! ! The following variables are currently hard-wired as snapshots ! (instantaneous rather than time-averages): -! divu, shear, sig1, sig2, trsig, mlt_onset, frz_onset, hisnap, aisnap +! divu, shear, sig1, sig2, sigI, trsig, mlt_onset, frz_onset, hisnap, aisnap ! ! Options for histfreq: '1','h','d','m','y','x', where x means that ! output stream will not be used (recommended for efficiency). @@ -227,6 +227,7 @@ module ice_history_shared f_strength = 'm', & f_divu = 'm', f_shear = 'm', & f_sig1 = 'm', f_sig2 = 'm', & + f_sigI = 'm', & f_dvidtt = 'm', f_dvidtd = 'm', & f_daidtt = 'm', f_daidtd = 'm', & f_dagedtt = 'm', f_dagedtd = 'm', & @@ -314,6 +315,7 @@ module ice_history_shared f_strength, & f_divu, f_shear , & f_sig1, f_sig2 , & + f_sigI, & f_dvidtt, f_dvidtd , & f_daidtt, f_daidtd , & f_dagedtt, f_dagedtd , & @@ -419,6 +421,7 @@ module ice_history_shared n_strength , & n_divu , n_shear , & n_sig1 , n_sig2 , & + n_sigI , & n_dvidtt , n_dvidtd , & n_daidtt , n_daidtd , & n_dagedtt , n_dagedtd , & diff --git a/cicecore/cicedynB/dynamics/ice_dyn_eap.F90 b/cicecore/cicedynB/dynamics/ice_dyn_eap.F90 index b2a528b9d..7510ce37e 100644 --- a/cicecore/cicedynB/dynamics/ice_dyn_eap.F90 +++ b/cicecore/cicedynB/dynamics/ice_dyn_eap.F90 @@ -93,7 +93,7 @@ subroutine eap (dt) cosw, sinw, denom1, uvel_init, vvel_init, arlx1i, & evp_prep1, evp_prep2, stepu, evp_finish, & basal_stress_coeff, basalstress - use ice_flux, only: rdg_conv, rdg_shear, prs_sig, strairxT, strairyT, & + use ice_flux, only: rdg_conv, rdg_shear, strairxT, strairyT, & strairx, strairy, uocn, vocn, ss_tltx, ss_tlty, iceumask, fm, & strtltx, strtlty, strocnx, strocny, strintx, strinty, & strocnxT, strocnyT, strax, stray, & @@ -176,7 +176,6 @@ subroutine eap (dt) rdg_shear(i,j,iblk) = c0 divu (i,j,iblk) = c0 shear(i,j,iblk) = c0 - prs_sig(i,j,iblk) = c0 e11(i,j,iblk) = c0 e12(i,j,iblk) = c0 e22(i,j,iblk) = c0 @@ -394,7 +393,6 @@ subroutine eap (dt) yieldstress11 (:,:,iblk), & yieldstress12 (:,:,iblk), & yieldstress22 (:,:,iblk), & - prs_sig (:,:,iblk), & rdg_conv (:,:,iblk), rdg_shear (:,:,iblk), & strtmp (:,:,:)) ! call ice_timer_stop(timer_tmp1) ! dynamics @@ -486,8 +484,8 @@ subroutine eap (dt) if ( basalstress ) then !$OMP PARALLEL DO PRIVATE(iblk) do iblk = 1, nblocks - taubx(:,:,iblk) = Cbu(:,:,iblk)*uvel(:,:,iblk) - tauby(:,:,iblk) = Cbu(:,:,iblk)*vvel(:,:,iblk) + taubx(:,:,iblk) = -Cbu(:,:,iblk)*uvel(:,:,iblk) + tauby(:,:,iblk) = -Cbu(:,:,iblk)*vvel(:,:,iblk) enddo !$OMP END PARALLEL DO endif @@ -1145,7 +1143,6 @@ subroutine stress_eap (nx_block, ny_block, & yieldstress11, & yieldstress12, & yieldstress22, & - prs_sig, & rdg_conv, rdg_shear, & strtmp) @@ -1196,7 +1193,6 @@ subroutine stress_eap (nx_block, ny_block, & real (kind=dbl_kind), dimension (nx_block,ny_block), & intent(inout) :: & - prs_sig , & ! replacement pressure, for stress calc shear , & ! strain rate II component (1/s) divu , & ! strain rate I component, velocity divergence (1/s) e11 , & ! components of strain rate tensor (1/s) @@ -1353,8 +1349,6 @@ subroutine stress_eap (nx_block, ny_block, & e22(i,j) = p5*p25*(divune + divunw + divuse + divusw - & tensionne - tensionnw - tensionse - tensionsw) * tarear(i,j) - prs_sig(i,j) = strength(i,j) - !----------------------------------------------------------------- ! elastic relaxation, see Eq. A12-A14 !----------------------------------------------------------------- diff --git a/cicecore/cicedynB/dynamics/ice_dyn_evp.F90 b/cicecore/cicedynB/dynamics/ice_dyn_evp.F90 index 18abfb04f..eb211681f 100644 --- a/cicecore/cicedynB/dynamics/ice_dyn_evp.F90 +++ b/cicecore/cicedynB/dynamics/ice_dyn_evp.F90 @@ -80,7 +80,7 @@ subroutine evp (dt) use ice_blocks, only: block, get_block, nx_block, ny_block use ice_domain, only: nblocks, blocks_ice, halo_info, maskhalo_dyn use ice_domain_size, only: max_blocks, ncat - use ice_flux, only: rdg_conv, rdg_shear, prs_sig, strairxT, strairyT, & + use ice_flux, only: rdg_conv, rdg_shear, strairxT, strairyT, & strairx, strairy, uocn, vocn, ss_tltx, ss_tlty, iceumask, fm, & strtltx, strtlty, strocnx, strocny, strintx, strinty, & strocnxT, strocnyT, strax, stray, & @@ -177,7 +177,6 @@ subroutine evp (dt) rdg_shear(i,j,iblk) = c0 divu (i,j,iblk) = c0 shear(i,j,iblk) = c0 - prs_sig(i,j,iblk) = c0 enddo enddo @@ -352,7 +351,6 @@ subroutine evp (dt) stress12_1(:,:,iblk), stress12_2(:,:,iblk), & stress12_3(:,:,iblk), stress12_4(:,:,iblk), & shear (:,:,iblk), divu (:,:,iblk), & - prs_sig (:,:,iblk), & rdg_conv (:,:,iblk), rdg_shear (:,:,iblk), & strtmp (:,:,:) ) ! endif ! yield_curve @@ -420,8 +418,8 @@ subroutine evp (dt) if ( basalstress ) then !$OMP PARALLEL DO PRIVATE(iblk) do iblk = 1, nblocks - taubx(:,:,iblk) = Cbu(:,:,iblk)*uvel(:,:,iblk) - tauby(:,:,iblk) = Cbu(:,:,iblk)*vvel(:,:,iblk) + taubx(:,:,iblk) = -Cbu(:,:,iblk)*uvel(:,:,iblk) + tauby(:,:,iblk) = -Cbu(:,:,iblk)*vvel(:,:,iblk) enddo !$OMP END PARALLEL DO endif @@ -551,7 +549,6 @@ subroutine stress (nx_block, ny_block, & stress12_1, stress12_2, & stress12_3, stress12_4, & shear, divu, & - prs_sig, & rdg_conv, rdg_shear, & str ) @@ -588,7 +585,6 @@ subroutine stress (nx_block, ny_block, & real (kind=dbl_kind), dimension (nx_block,ny_block), & intent(inout) :: & - prs_sig , & ! replacement pressure, for stress calc shear , & ! strain rate II component (1/s) divu , & ! strain rate I component, velocity divergence (1/s) rdg_conv , & ! convergence term for ridging (1/s) @@ -698,7 +694,6 @@ subroutine stress (nx_block, ny_block, & c0nw = strength(i,j)/max(Deltanw,tinyarea(i,j)) c0sw = strength(i,j)/max(Deltasw,tinyarea(i,j)) c0se = strength(i,j)/max(Deltase,tinyarea(i,j)) - prs_sig(i,j) = c0ne*Deltane ! northeast c1ne = c0ne*arlx1i c1nw = c0nw*arlx1i diff --git a/cicecore/cicedynB/dynamics/ice_dyn_shared.F90 b/cicecore/cicedynB/dynamics/ice_dyn_shared.F90 index 44185ec8c..71729e1d0 100644 --- a/cicecore/cicedynB/dynamics/ice_dyn_shared.F90 +++ b/cicecore/cicedynB/dynamics/ice_dyn_shared.F90 @@ -851,6 +851,8 @@ subroutine basal_stress_coeff (nx_block, ny_block, icellu, & uold, vold, & Cbu) + use ice_constants, only: c0, c1 + integer (kind=int_kind), intent(in) :: & nx_block, ny_block, & ! block dimensions icellu ! no. of cells where icetmask = 1 @@ -878,7 +880,7 @@ subroutine basal_stress_coeff (nx_block, ny_block, icellu, & hu, & ! volume per unit area of ice at u location (mean thickness) hwu, & ! water depth at u location hcu, & ! critical thickness at u location - k1 = 8.0_dbl_kind , & ! first free parameter for landfast parametrization + k1 = 80.0_dbl_kind , & ! first free parameter for landfast parametrization k2 = 15.0_dbl_kind, & ! second free parameter (Nm^-3) for landfast parametrization u0 = 5e-5_dbl_kind, & ! residual velocity (m/s) CC = 20.0_dbl_kind ! CC=Cb factor in Lemieux et al 2015 @@ -900,11 +902,9 @@ subroutine basal_stress_coeff (nx_block, ny_block, icellu, & hcu = au * hwu / k1 ! 2- calculate stress factor - if (au > p01 .and. hu > hcu ) then - ! endif - Cbu(i,j) = ( k2 / (sqrt(uold(i,j)**2 + vold(i,j)**2) + u0) ) & - * (hu - hcu) * exp(-CC * (1 - au)) - endif + + Cbu(i,j) = ( k2 / (sqrt(uold(i,j)**2 + vold(i,j)**2) + u0) ) & + * max(c0,(hu - hcu)) * exp(-CC * (c1 - au)) enddo ! ij @@ -919,8 +919,9 @@ end subroutine basal_stress_coeff subroutine principal_stress(nx_block, ny_block, & stressp_1, stressm_1, & - stress12_1, prs_sig, & - sig1, sig2) + stress12_1, strength, & + sig1, sig2, & + sigI) use ice_constants, only: spval_dbl, p5, c4 @@ -931,11 +932,12 @@ subroutine principal_stress(nx_block, ny_block, & stressp_1 , & ! sigma11 + sigma22 stressm_1 , & ! sigma11 - sigma22 stress12_1, & ! sigma12 - prs_sig ! replacement pressure, for stress calc + strength ! for normalization of sig1 and sig2 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(out):: & - sig1 , & ! principal stress component - sig2 ! principal stress component + sig1 , & ! normalized principal stress component + sig2 , & ! normalized principal stress component + sigI ! internal ice pressure (N/m) ! local variables @@ -950,16 +952,20 @@ subroutine principal_stress(nx_block, ny_block, & do j = 1, ny_block do i = 1, nx_block - if (prs_sig(i,j) > puny) then + if (strength(i,j) > puny) then + ! not normalized yet sig1(i,j) = (p5*(stressp_1(i,j) & - + sqrt(stressm_1(i,j)**2+c4*stress12_1(i,j)**2))) & - / prs_sig(i,j) + + sqrt(stressm_1(i,j)**2+c4*stress12_1(i,j)**2))) sig2(i,j) = (p5*(stressp_1(i,j) & - - sqrt(stressm_1(i,j)**2+c4*stress12_1(i,j)**2))) & - / prs_sig(i,j) + - sqrt(stressm_1(i,j)**2+c4*stress12_1(i,j)**2))) + sigI(i,j) = -p5*( sig1(i,j) + sig2(i,j) ) + ! normalization of sig1 and sig2 + sig1(i,j) = sig1(i,j) / strength(i,j) + sig2(i,j) = sig2(i,j) / strength(i,j) else sig1(i,j) = spval_dbl sig2(i,j) = spval_dbl + sigI(i,j) = spval_dbl endif enddo enddo diff --git a/cicecore/cicedynB/general/ice_flux.F90 b/cicecore/cicedynB/general/ice_flux.F90 index 15324491c..fce33c19f 100644 --- a/cicecore/cicedynB/general/ice_flux.F90 +++ b/cicecore/cicedynB/general/ice_flux.F90 @@ -58,8 +58,9 @@ module ice_flux ! diagnostic real (kind=dbl_kind), dimension (nx_block,ny_block,max_blocks), public :: & - sig1 , & ! principal stress component - sig2 , & ! principal stress component + sig1 , & ! normalized principal stress component + sig2 , & ! normalized principal stress component + sigI , & ! internal ice pressure (N/m) taubx , & ! basal stress (x) (N/m^2) tauby , & ! basal stress (y) (N/m^2) strairx , & ! stress on ice by air, x-direction @@ -108,7 +109,6 @@ module ice_flux ! internal real (kind=dbl_kind), dimension (nx_block,ny_block,max_blocks), public :: & - prs_sig , & ! replacement pressure, for stress calc fm , & ! Coriolis param. * mass in U-cell (kg/s) Cbu ! coefficient for basal stress (landfast ice) @@ -721,7 +721,6 @@ subroutine init_history_dyn if (tr_iage) & dagedtd (:,:,:) = trcr(:,:,nt_iage,:) ! temporary initial age fm (:,:,:) = c0 - prs_sig (:,:,:) = c0 ardgn (:,:,:,:) = c0 vrdgn (:,:,:,:) = c0 krdgn (:,:,:,:) = c1 diff --git a/cicecore/cicedynB/infrastructure/io/io_binary/ice_history_write.F90 b/cicecore/cicedynB/infrastructure/io/io_binary/ice_history_write.F90 index 863d7314c..2fc173bae 100644 --- a/cicecore/cicedynB/infrastructure/io/io_binary/ice_history_write.F90 +++ b/cicecore/cicedynB/infrastructure/io/io_binary/ice_history_write.F90 @@ -147,7 +147,7 @@ subroutine ice_write_hist(ns) if (histfreq(ns) == '1' .or. .not. hist_avg & .or. n==n_divu(ns) .or. n==n_shear(ns) & ! snapshots .or. n==n_sig1(ns) .or. n==n_sig2(ns) & - .or. n==n_trsig(ns) & + .or. n==n_sigI(ns) .or. n==n_trsig(ns) & .or. n==n_mlt_onset(ns) .or. n==n_frz_onset(ns) & .or. n==n_hisnap(ns) .or. n==n_aisnap(ns)) then write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), & diff --git a/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_history_write.F90 b/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_history_write.F90 index 899302b6f..51ccd48dd 100644 --- a/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_history_write.F90 +++ b/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_history_write.F90 @@ -505,7 +505,8 @@ subroutine ice_write_hist (ns) !----------------------------------------------------------------- if (hist_avg) then if (TRIM(avail_hist_fields(n)%vname)/='sig1' & - .or.TRIM(avail_hist_fields(n)%vname)/='sig2') then + .or.TRIM(avail_hist_fields(n)%vname)/='sig2' & + .or.TRIM(avail_hist_fields(n)%vname)/='sigI') then status = nf90_put_att(ncid,varid,'cell_methods','time: mean') if (status /= nf90_noerr) call abort_ice( & 'Error defining cell methods for '//avail_hist_fields(n)%vname) @@ -515,7 +516,7 @@ subroutine ice_write_hist (ns) if (histfreq(ns) == '1' .or. .not. hist_avg & .or. n==n_divu(ns) .or. n==n_shear(ns) & ! snapshots .or. n==n_sig1(ns) .or. n==n_sig2(ns) & - .or. n==n_trsig(ns) & + .or. n==n_sigI(ns) .or. n==n_trsig(ns) & .or. n==n_mlt_onset(ns) .or. n==n_frz_onset(ns) & .or. n==n_hisnap(ns) .or. n==n_aisnap(ns)) then status = nf90_put_att(ncid,varid,'time_rep','instantaneous') diff --git a/cicecore/cicedynB/infrastructure/io/io_pio/ice_history_write.F90 b/cicecore/cicedynB/infrastructure/io/io_pio/ice_history_write.F90 index 5c2e8b7d0..d87c306b6 100644 --- a/cicecore/cicedynB/infrastructure/io/io_pio/ice_history_write.F90 +++ b/cicecore/cicedynB/infrastructure/io/io_pio/ice_history_write.F90 @@ -425,7 +425,8 @@ subroutine ice_write_hist (ns) ! Add cell_methods attribute to variables if averaged if (hist_avg .and. histfreq(ns) /= '1') then if (TRIM(avail_hist_fields(n)%vname)/='sig1' & - .or.TRIM(avail_hist_fields(n)%vname)/='sig2') then + .or.TRIM(avail_hist_fields(n)%vname)/='sig2' & + .or.TRIM(avail_hist_fields(n)%vname)/='sigI') then status = pio_put_att(File,varid,'cell_methods','time: mean') endif endif @@ -433,7 +434,7 @@ subroutine ice_write_hist (ns) if (histfreq(ns) == '1' .or. .not. hist_avg & .or. n==n_divu(ns) .or. n==n_shear(ns) & ! snapshots .or. n==n_sig1(ns) .or. n==n_sig2(ns) & - .or. n==n_trsig(ns) & + .or. n==n_sigI(ns) .or. n==n_trsig(ns) & .or. n==n_mlt_onset(ns) .or. n==n_frz_onset(ns) & .or. n==n_hisnap(ns) .or. n==n_aisnap(ns)) then status = pio_put_att(File,varid,'time_rep','instantaneous') diff --git a/configuration/scripts/ice_in b/configuration/scripts/ice_in index f543b2e57..08a5b082f 100644 --- a/configuration/scripts/ice_in +++ b/configuration/scripts/ice_in @@ -420,6 +420,7 @@ f_shear = 'm' f_sig1 = 'm' f_sig2 = 'm' + f_sigI = 'm' f_dvidtt = 'm' f_dvidtd = 'm' f_daidtt = 'm'