source: issm/trunk/src/m/classes/lovenumbers.py

Last change on this file was 28013, checked in by Mathieu Morlighem, 16 months ago

merged trunk-jpl and trunk for revision 28011
File size: 8.4 KB
Line 
1import numpy as np
2from checkfield import checkfield
3from fielddisplay import fielddisplay
4from getlovenumbers import getlovenumbers
5from pairoptions import pairoptions
6from WriteData import WriteData
7
8
9class lovenumbers(object): #{{{
10 """lovenumbers class definition
11
12 Usage:
13 lovenumbers = lovenumbers()
14 lovenumbers = lovenumbers('maxdeg', 10000, 'referenceframe', 'CF');
15
16 Choose numbers of degrees required (1000 by default) and reference frame
17 (between CF and CM; CM by default)
18 """
19
20 def __init__(self, *args): #{{{
21 # Loading love numbers
22 self.h = [] # Provided by PREM model
23 self.k = [] # idem
24 self.l = [] # idem
25
26 # Tidal love numbers for computing rotational feedback
27 self.th = []
28 self.tk = []
29 self.tl = []
30 self.tk2secular = 0 # deg 2 secular number
31 self.pmtf_colinear = []
32 self.pmtf_ortho = []
33
34 # Time/frequency for visco-elastic love numbers
35 self.timefreq = []
36 self.istime = 1
37
38 options = pairoptions(*args)
39 maxdeg = options.getfieldvalue('maxdeg', 1000)
40 referenceframe = options.getfieldvalue('referenceframe', 'CM')
41 self.setdefaultparameters(maxdeg, referenceframe)
42 # }}}
43
44 def __repr__(self): #{{{
45 s = ' lovenumbers parameters:\n'
46 s += '{}\n'.format(fielddisplay(self, 'h', 'load Love number for radial displacement'))
47 s += '{}\n'.format(fielddisplay(self, 'k', 'load Love number for gravitational potential perturbation'))
48 s += '{}\n'.format(fielddisplay(self, 'l', 'load Love number for horizontal displacements'))
49 s += '{}\n'.format(fielddisplay(self, 'th', 'tidal load Love number (deg 2)'))
50 s += '{}\n'.format(fielddisplay(self, 'tk', 'tidal load Love number (deg 2)'))
51 s += '{}\n'.format(fielddisplay(self, 'tl', 'tidal load Love number (deg 2)'))
52 s += '{}\n'.format(fielddisplay(self, 'tk2secular', 'secular fluid Love number'))
53 s += '{}\n'.format(fielddisplay(self, 'pmtf_colinear', 'Colinear component of the Polar Motion Transfer Function (e.g. x-motion due to x-component perturbation of the inertia tensor)'))
54 s += '{}\n'.format(fielddisplay(self, 'pmtf_ortho', 'Orthogonal component of the Polar Motion Transfer Function (couples x and y components, only used for Chandler Wobble)'))
55 s += '{}\n'.format(fielddisplay(self, 'istime', 'time (default: 1) or frequency love numbers (0)'))
56 s += '{}\n'.format(fielddisplay(self, 'timefreq', 'time/frequency vector (yr or 1/yr)'))
57 s += '{}\n'.format(fielddisplay(self, 'pmtf_colinear', 'Colinear component of the Polar Motion Transfer Function (e.g. x-motion due to x-component perturbation of the inertia tensor)'))
58 s += '{}\n'.format(fielddisplay(self, 'pmtf_ortho', 'Orthogonal component of the Polar Motion Transfer Function (couples x and y components, only used for Chandler Wobble)'))
59 return s
60 # }}}
61
62 def setdefaultparameters(self, maxdeg, referenceframe): #{{{
63 # Initialize love numbers
64 self.h = getlovenumbers('type', 'loadingverticaldisplacement', 'referenceframe', referenceframe, 'maxdeg', maxdeg).reshape(-1,1)
65 self.k = getlovenumbers('type', 'loadinggravitationalpotential', 'referenceframe', referenceframe, 'maxdeg', maxdeg).reshape(-1,1)
66 self.l = getlovenumbers('type', 'loadinghorizontaldisplacement', 'referenceframe', referenceframe, 'maxdeg', maxdeg).reshape(-1,1)
67 self.th = getlovenumbers('type', 'tidalverticaldisplacement', 'referenceframe', referenceframe, 'maxdeg', maxdeg).reshape(-1,1)
68 self.tk = getlovenumbers('type', 'tidalgravitationalpotential', 'referenceframe', referenceframe, 'maxdeg', maxdeg).reshape(-1,1)
69 self.tl = getlovenumbers('type', 'tidalhorizontaldisplacement', 'referenceframe', referenceframe, 'maxdeg', maxdeg).reshape(-1,1)
70
71 # Secular fluid love number
72 self.tk2secular = 0.942
73
74 self.pmtf_colinear = np.array([0.0]).reshape(-1, 1)
75 self.pmtf_ortho = np.array([0.0]).reshape(-1, 1)
76 if maxdeg >= 2:
77 self.pmtf_colinear = ((1.0 + self.k[2, :]) / (1.0 - self.tk[2, :] / self.tk2secular)).reshape(-1, 1) # Valid only for elastic regime, not viscous. Also neglects chandler wobble.
78 self.pmtf_ortho = np.array([0.0]).reshape(-1, 1)
79 # Time
80 self.istime = 1 # Temporal love numbers by default
81 self.timefreq = np.zeros(1) # Elastic case by default
82 return self
83 # }}}
84
85 def checkconsistency(self, md, solution, analyses): #{{{
86 if ('SealevelchangeAnalysis' not in analyses) or (solution == 'TransientSolution' and not md.transient.isslc):
87 return
88 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.h', 'NaN', 1, 'Inf', 1)
89 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.k', 'NaN', 1, 'Inf', 1)
90 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.l', 'NaN', 1, 'Inf', 1)
91
92 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.th', 'NaN', 1, 'Inf', 1)
93 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.tk', 'NaN', 1, 'Inf', 1)
94 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.tl', 'NaN', 1, 'Inf', 1)
95 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.tk2secular', 'NaN', 1, 'Inf', 1)
96 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.pmtf_colinear', 'NaN', 1, 'Inf', 1)
97 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.pmtf_ortho', 'NaN', 1, 'Inf', 1)
98 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.timefreq', 'NaN', 1, 'Inf', 1)
99 md = checkfield(md, 'fieldname', 'solidearth.lovenumbers.istime', 'NaN', 1, 'Inf', 1, 'values', [0, 1])
100
101 # Check that love numbers are provided at the same level of accuracy
102 if (self.h.shape[0] != self.k.shape[0]) or (self.h.shape[0] != self.l.shape[0]):
103 raise ValueError('lovenumbers error message: love numbers should be provided at the same level of accuracy')
104
105 ntf = len(self.timefreq)
106 if (np.shape(self.h)[1] != ntf or np.shape(self.k)[1] != ntf or np.shape(self.l)[1] != ntf or np.shape(self.th)[1] != ntf or np.shape(self.tk)[1] != ntf or np.shape(self.tl)[1] != ntf or np.shape(self.pmtf_colinear)[1] != ntf or np.shape(self.pmtf_ortho)[1] != ntf):
107 raise ValueError('lovenumbers error message: love numbers should have as many time/frequency steps as the time/frequency vector')
108
109 if self.istime and self.timefreq[0] != 0:
110 raise ValueError('temporal love numbers must start with elastic response, i.e. timefreq[0] = 0')
111 return md
112 # }}}
113
114 def defaultoutputs(self, md): #{{{
115 return[]
116 # }}}
117
118 def marshall(self, prefix, md, fid): #{{{
119 WriteData(fid, prefix, 'object', self, 'fieldname', 'h', 'name', 'md.solidearth.lovenumbers.h', 'format', 'DoubleMat', 'mattype', 1)
120 WriteData(fid, prefix, 'object', self, 'fieldname', 'k', 'name', 'md.solidearth.lovenumbers.k', 'format', 'DoubleMat', 'mattype', 1)
121 WriteData(fid, prefix, 'object', self, 'fieldname', 'l', 'name', 'md.solidearth.lovenumbers.l', 'format', 'DoubleMat', 'mattype', 1)
122
123 WriteData(fid, prefix, 'object', self, 'fieldname', 'th', 'name', 'md.solidearth.lovenumbers.th', 'format', 'DoubleMat', 'mattype', 1)
124 WriteData(fid, prefix, 'object', self, 'fieldname', 'tk', 'name', 'md.solidearth.lovenumbers.tk', 'format', 'DoubleMat', 'mattype', 1)
125 WriteData(fid, prefix, 'object', self, 'fieldname', 'tl', 'name', 'md.solidearth.lovenumbers.tl', 'format', 'DoubleMat', 'mattype', 1)
126 WriteData(fid, prefix, 'object', self, 'fieldname', 'pmtf_colinear','name','md.solidearth.lovenumbers.pmtf_colinear','format','DoubleMat','mattype',1)
127 WriteData(fid, prefix, 'object', self, 'fieldname', 'pmtf_ortho','name','md.solidearth.lovenumbers.pmtf_ortho','format','DoubleMat','mattype',1)
128 WriteData(fid, prefix, 'object', self, 'data', self.tk2secular, 'fieldname', 'lovenumbers.tk2secular', 'format', 'Double')
129
130 if (self.istime):
131 scale = md.constants.yts
132 else:
133 scale = 1.0 / md.constants.yts
134 WriteData(fid, prefix, 'object', self, 'fieldname', 'istime', 'name', 'md.solidearth.lovenumbers.istime', 'format', 'Boolean')
135 WriteData(fid, prefix, 'object', self, 'fieldname', 'timefreq', 'name', 'md.solidearth.lovenumbers.timefreq', 'format', 'DoubleMat', 'mattype', 1, 'scale', scale);
136 # }}}
137
138 def extrude(self, md): #{{{
139 return
140 # }}}
Note: See TracBrowser for help on using the repository browser.