Fortran programs

gotsalsm / gotsalzm


USAGE:
gotsal{s|z}m @ <ins-file> [:<target>]

PURPOSE:
Compute global loading field; also, iterate self-attraction and -loading of ocean tide generating potential
Uses netCDF tide fields for input and output but in a non-orthodox way (header information not consistent with what GMT-grd programs expect).
 
Equilibrium ocean response may be computed (using subroutine SET_ZM_SH, oteu0162.f). Note that the default settings imply unit-normalised
spherical harmonics, where the P20 term is sign-reversed. The Doodson/Bartels normalizsation to obtain a positive maximum of unity is optional.

gotsalzm uses complex-valued Greens functions and offers many more options (marked Z in Namelist table below)
gotsalzm employs Chebychev-Gauss integration while gotsalsm uses brute-force summation of regular sub-cells.

Chebychev-Gauss may consume significant computation time; careful trading-off of the criteria for the iteration of the integral is encouraged. On the other hand, previous use of the brute-force method did not anticipate the comparably large errors incurred by to coarse sub-boxing.

FUNCTION:
Accelerated convolution, integrated loading kernel in near zone.

SOURCE CODE:
Oload/afor/p/m/gotsalsm.f
Oload/afor/p/mz/gotsalzm.f

INSTRUCTION FILE:
    Namelist &param
    Open-file block

    Namelist parameters
Parameter
- type -
Explanation
Default


comp



- char*4  -

Input / excitation
.
Component of the loading effect to be computed,
any of OTEP, RADI, GRAV, PTAN

Special value $ENV : obtain comp from environment GOTSAL_COMP



OTEP

gotd_path
- char*64 -
Path to the input netCDF directory
/home/loading/tide-models/netCDF/
gotd_name
- char*64 -
Tide model sub-directory
Special value $ENV : obtain gotd_name from environment GOTSAL_MODEL
FES94.1
tide
- char*4  -
Darwin symbol
The tide file names are composed as follows:
gotd_path//gotd_name//'/amp.'//tide
gotd_path//gotd_name//'/pha.'//tide
Special value $ENV : obtain tide from environment GOTSAL_TIDE

M2
qeqlt
- logical -
Use the land-sea distribution of the input tide model,
fill the array at its wet nodes with a spherical harmonic
.false.
ndeg_eqlt
mord_eqlt
lcs_eqlt
Z
- integer -
With qeqlt=.true., spherical harmonic degree and ...
order of excitation of equilibrium tide

Z and a selector 0 to 3 for cos, i sin, complex exp, real sin 
2
0
Z0
zmass
mdistrib_flag
q_from_zero
- complex -
- integer -
- logical -
Distribute a mass [kg] at flag value mdistrib_flag
(requires input of a flag array from unit 21 2)
wipe the input field before applying the mass.
(Currently, the only option is a uniform slab.)

(0.0d0,0.0d0)
2
.false.


dir_out



- char*64 -

Output
.
Path to the output netCDF files.
The files will be named
dir_out//comp//'_amp.'//tide
dir_out//comp//'_pha.'//tide
Special value $ENV : obtain dir_out from environment GOTSAL_DIROUT


./

tide_out
- char*4  -
Z If non-empty, the extension of the output file name.
Default is the name of the input tide.

' '
qout_add_extpot
- logical - Z Under OTEP and qeqlt=.true., add the external potential to the loading effect.
.false.
qextpot_out
- logical - Z Under qeqlt=.true., output the external potential field. The file names
will be marked with 'EXTF'
.false.
qunity_scale_out
- logical - Z Under qeqlt=.true., scale the output relative to the maximum (not max abs!)
of the forcing field.

.false.
qout_add_mass
- logical - (not Z) Under OTEP and |zmass| > 1.0 kg, add the loading masses to the OTEP field at wet points and impose them on land. That's to lay a finishing hand to an iteration of OTEP. If you select this option, prefer to set
qout_allnodes=.false. and eventually combine the nodes that carry the
mdistrib_flag with the wet nodes:
nchange_flags=1, change_flags(1,1)=2, change_flags(1,2)=1

 .true.
qout_allnodes


nchange_flags
change_flags(10,2)
- logical -


- integer -
- integer -
.true. is the preferred option for global maps of loading effects.
.false. in conjunction with iteration of OTEP and producing self-consistent maps of sea level.
Re-flagging option for output. For instance, in the case of |zmass| > 1 kg, re-flag nodes change_flag(j,1) as
change_flag(j,2)
for 1 <= j <= nchange_flags
.true.


0
20*0
qversal
- logical -
In output file name (tide), use all capitals in the Darwin symbols of diurnal, semi-durnal etc. tides

.false.
kcentered
- integer -
Output grid can be recentered on Europe or Pacific; or left as it was.
0 - don't change
1 - Europe
2 - Pacific
0
long_widen
- integer -
On output, add this number of meridional stripes to the east and west,
creating a cylindrical overlap. (The subroutines are in extend-ew.f)

0
kind_of_output(2)

- char*4 -
'C...' (e.g. 'CMPX') will produce files _real.tide and _imag.tide
'P...' (e.g. 'POLR') will produce files _amp.tide and _pha.tide
'....' will use the default (or current setting) in putnetcdf (getnetcdf.f)
'*...' will skip output.
_amp & _pha.tide


'.','*'


fn_greens



- char*32 -

Greens function

.
Greens function name


mc00egbo

dir_greens
- char*32 - Greens function directory
The file name is composed ass follows:
dir_greens//fn_greens//'.dat'
Greenstb/
lfgo
lfso, lfsi
- integer -
If > 0 and file open, list the normalized expanded Greens function.
Write integrated Greens values to lfso, (re-)input them on lfsi.
41
42,43
x0
- real*8 -
Lower distance limit in degrees for Greens function expansion
1)       0.01
ng
- integer -
Length of array for Greens function expansion

1)       3000
qgreen_real



qgotsa
l_real
- logical -



- logical -

Z Take only the real-part of the Greens function (a parameter to expglgz).
   The Greens function listed on file (unit lfgo) has imaginary part
   zero as will the stored integrated values on file (unit lfso).

  Subroutine gotsal (gotsalzs.f) uses only the real-part of Greens in the
   convolution. Complex values are written to lfgo and lfso
  
         .false.




.false.


q_iterate


- logical -
Self-attraction and -loading iteration
.
Iterate the OTEP solution


.false.
n_iterate
- integer -
Maximum number of iterations
1
crit_iterate
- real*8  -
Continue iteration while the updates differ less than this value
(relative measure)

0.01
qout_add_extpot
qout_add_mass

- logical -
cf Output section above
.false.


crtdis
qclose_by_nr

qreciproc_acc 


- real*8  -
- logical -
- logical -
Convolution and loading-kernel integration
.
Maximum distance in degrees for integration of Greens function
Z Instead of crtdis citerion, integrate only the eight neighbours
Z Apply increased accuracy in ring-to-ring reciprocity

  

  
1
)      10.0
    .false.
  .false.
eps
eps_s _i _n 
- real*8  -
- real*8  -
(not Z) Convergence criterion for integration
Z      ... cases self, not self, sloppy

1)     0.001
 
3*0.001
nresi
- integer - Number of sub-intervals in integration of neighbour cells
1)         9
nress
- integer - Starting number of sub-intervals in integration of neighbour cells 1)         3
maxstp.
maxstp_s _i _n
- integer -
- integer -
(not Z) Break iteration when number of sub-intervals exceeds this number
Z Break iteration, cases self, not self, sloppy

.4500
 
8200,1100,1100


ipgre


- integer -
Miscellaneous

Print out Greens table at this interval


10
mod_print
- integer - In the gotsal subroutine, diagnostic print:
row number modulo this number.
mod_print < 90000: print the account of convolution (how many source rings contributed to every target ring). Also, report of the input rings' RMS at the stage of FFT.
99999
check_at_flags(10)
- integer -
In the context of loading at specific types of flags (cf zmass), print a statistic for each value of this array. A value -1 terminates the array.

10*-1
qprint_bsmean
- logical - Z Monitor Greens integration in self case
 
.true.
greensp_monlat
iun_greensp_mon
- real*8  -
- integer -
Z Output the Greens function spectrum at this latitude to ...
  ... this file unit (binary output with subroutine outsp)

100.0 i.e. don't
-1
qno_self
qbox_on_box
- logical - Z bsmean is skipped in self-case.
  if .false., bsmean is executed in self-case but result set to zero

  .false.
 .true.
sloppy_above
- real*8  - Z Relax the iteration criterion above this latitude

100.0
gint_monitor_lat
- real*8  - Z Target latitude to monitor Greens integration in non-self cases

100.0
ntest_bsmean
- integer -
Z Test the box integration of Greens with this value as the starting
   number of nodes. gotsal will stop after the test.
 -1
1) Defaults are assigned in Block data program bgotsa.f . All other are assigned in Main Data statement.
2) Program loadland produces such an array (Oload/afor/p/m/loadland.f)

    Open-file block
Unit
   format
purpose
21
binary
Optional; a flag array required for solving problems of loading masses in specific areas
41
 ASCII 
Greens function tabulated at interval ipgre
42
 binary 
If open, input an existing table of integrated Greens function; must pertain to the component comp
43
binary
If 42 not open, output the integrated Greens function

  EXAMPLES:

  1)  ins-file:

IT>
 &param

 tide='M2'
 eps=2.d-4, nress=7, maxstp=4500
 fn_greens='mc00egbcw'
 NG=30000
 X0=0.001
 dir_out='FES94/', qout_allnodes=.true.
 q_iterate=.true., n_iterate=2, crit_iterate=0.001
 &end
41 B Greenstb/hdg/expgs.otep
42 ^ Greenstb/hdg/btbls.otep
43 * < Greenstb/hdg/btbls.otep
   Q

Dwells on defaults for comp='OTEP' and gotd_name='FES94.1'

Exec:

    gotsalsm @ gotsalsm.ins :IT

Reads a previously created file with integrated Greens kernel values (unit 42). Subdirectory hdg collects half-degree items (qdg - quarter-degree, 8dg - eight nodes per degree (future))

  2)  Mass loading problem: Greenland ice discharge
        Prepare:
mkdir -p greenland
plot-gmtcoast g
#
There's an easier method using GoogleEarth -> draw path -> save-as kml
# kml2lola
<kml-file>
#
take care of the reported xy-file
#
mv gmt-segments.xy tmp/greenland.xy
loadland @ loadland.ins
cd greenland
ln -s amp.Gi _amp.gi
ln -s pha.Gi _pha.gi
mkdir -p unmasked with-load

# File gotsalsm-qdg-otep.ins :

setenv GOTSAL_TIDE   GI

setenv GOTSAL_MODEL  greenland/_
setenv GOTSAL_DIROUT greenland
setenv GOTSAL_PATH ./
echo " "
exit

O>
 &param
 tide='$ENV'
 comp='OTEP'
 gotd_name='$ENV'
 gotd_path='$ENV'
 eps=2.d-4, nress=7, maxstp=4500
 fn_greens='mc00egbcw'
 NG=30000
 X0=0.001
 dir_out='$ENV', qout_allnodes=.false.
 change_flags(1,1)=2, change_flags(1,2)=1
 long_widen=0
 mod_print=10
 kcentered=0
 kstop=-1
 zmass=(250.d12,0.0d0)
 q_from_zero=.true.
 kind_of_output='CMPX'
 check_at_flags=0,1,2,-9
 q_iterate=.true., n_iterate=2
 &end
21 ^ greenland/flz.dat
41 B Greenstb/qdg/expgs.otep
42 ^ Greenstb/qdg/btbls.otep
43 * < Greenstb/qdg/btbls.otep
   Q

Exec:
source gotsalsm-qdg-otep.ins
gotsalsm @ gotsalsm-qdg-otep.ins :O | tee greenland/gotsalsm-it2.log
    Run a second time with  kind_of_output=' ' to produce amp and pha for loading calculation.
     For the different melting regions, associate different, faked  tide. Log in at loading, cd bin, and run
h[loading]~/bin> olfg @ olfg-grav.ins
cd greenland
ln -s OTEP_amp.Gi amp.m2
ln -s OTEP_pha.Gi pha.m2

grdmath OTEP_real.Gi -10  10 INRANGE 0 NAN OTEP_real.Gi MUL -1.0 MUL = with-load/OTEP_real.Gi
grdmath OTEP_amp.Gi  -10  10 INRANGE 0 NAN OTEP_amp.Gi  MUL = unmasked/OTEP_amp.Gi
grdmath OTEP_pha.Gi -400 400 INRANGE 0 NAN OTEP_pha.Gi  MUL = unmasked/OTEP_pha.Gi

# grdmath greenland/OTEP_real.Gi greenland/OTEP_imag.Gi HYPOT = greenland/OTEP_amp.Gi
plot-global-map -euro -Z -0.005,0.005 -mrk greenland greenland/with-load/OTEP_real.Gi

gotsalzm:
setenv GOTSAL_TIDE   M2
setenv GOTSAL_MODEL  GOT4.8/
setenv GOTSAL_DIROUT GOT4.8q
setenv GOTSAL_PATH '/home/loading/tide-models/netCDF/'
mkdir -p $GOTSAL_DIROUT
gotsalzm @ gotsalzm-prem-hdg-rt.ins :PR

PR>
 &param
 tide='$ENV'
 comp='RADI'
 gotd_name='$ENV'
 gotd_path='$ENV'
 eps=2.d-4, nress=128, nresi=128
 crtdis=7.
 fn_greens='mc24qprc'
 lfsi=42, lfso=43
 NG=300000
 X0=0.00001
 qallow_overrange=.true.
 dir_out='$ENV', qout_allnodes=.true.
 long_widen=0
 mod_print_conv=20, mod_print_fft=99999, qprint_bsmean=.true.
 kcentered=1
 maxstp_n=512, sloppy_above=80.0
 nresi=64, maxstp_i=256
 &end
41 B Greenstb/hdg/expgs.radi
42 * ^ Greenstb/hdg/btbls-premq.radi
43 < Greenstb/hdg/btbls-premq.radi
   Q


Equilibrium lunar nodal tide with a complex, anelastic Greens function:
source gotsalzm-8dg-otep.ins
gotsalzm @ gotsalzm-8dg-otep.ins :O > ! HAMTIDE$GRE$VAR/gotsalzm-8dg-otep-p20.log &
# for $ENV in gotd_name:
setenv VAR -s3deg
setenv GRE -PREMN15
setenv GOTSAL_MODEL  HAMTIDE/
setenv GOTSAL_DIROUT HAMTIDE$GRE$VAR
setenv GOTSAL_TIDE M2
setenv GOTSAL_PATH '/home/loading/tide-models/netCDF/'
mkdir -p $GOTSAL_DIROUT
exit
O>
 &param
 qeqlt=.true.,ndeg_eqlt=2,mord_eqlt=0
 lcs_eqlt=0,  tide_out='P20', qunity_scale_out=.true.
 tide='$ENV'
 comp='OTEP'
 gotd_name='$ENV'
 gotd_path='$ENV'
 eps_i=1.d-3, eps_s=1.e-3, nress=128, nresi=8
 maxstp_s=1200
 fn_greens='PREMN15'
 NG=30000
 X0=0.001
 dir_out='$ENV'
 qout_allnodes=.false.
 long_widen=0
 mod_print=10
 kcentered=0
 qextpot_out=.false.
 kstop=-1
 qclose_by_nr=.false., crtdis=3.0
 qgreen_real=.false., qgotsal_real=.false.
 q_iterate=.true., n_iterate=2
 qout_add_extpot=.true.
 &end
21 ^ greenland/flz.dat
41 B Greenstb/8dg/expgs_nodal.otep
42 * ^ Greenstb/8dg/btbls_nodal.otep
43 < Greenstb/8dg/btbls_nodal.otep
   Q