← Back to branch summary

~amcg-stokes/fluidity/multimaterial_diagnostic_dependencies

« back to all changes in this revision

Viewing changes to tests/viscous_fs_zhong_spatial_explicit_varrho/solution.py

  • Committer: Cian Wilson
  • Date: 2012-12-10 20:21:07 UTC
  • mfrom: (4132.1.7 fluidity)
  • Revision ID: cwilson@ldeo.columbia.edu-20121210202107-5wppwqcba4bfd1r3
Merging in changes from lp:fluidity.

Show diffs side-by-side

added added

removed removed

Lines of Context:
tests/viscous_fs_zhong_spatial_explicit_varrho/solution.py
 
1
 
 
2
from fluidity_tools import stat_parser as stat
 
3
from math import log, pi, cos, sin, cosh, sinh, exp
 
4
from numpy import sqrt, interp
 
5
from scipy.integrate import quad
 
6
 
 
7
wavelengthfactor = 1.0
 
8
 
 
9
eta0 = 0. # m
 
10
xi0 = 0. # m
 
11
#eta0 = 1.e3 # m
 
12
#xi0 = 1.e3 # m
 
13
depth=1.5e6 # m
 
14
rho0 = 4500. # kg/m**3
 
15
rhou = 2*rho0 # kg/m**3
 
16
mu0 = 1.e21  # Pas
 
17
g = 10.      # m/s**2
 
18
D = 3.e6     # m
 
19
deltaT = 2000 #K
 
20
alpha = 2.e-5 #/K
 
21
 
 
22
#eta0 = 1./3000. # m
 
23
#xi0 = 1./3000. # m
 
24
#depth=0.5 # m
 
25
#rho0 = 1.0 # kg/m**3
 
26
#rhou = 2.0 # kg/m**3
 
27
#mu0 = 1.0  # Pas
 
28
#g = 1.0      # m/s**2
 
29
#D = 1.0     # m
 
30
#deltaT = 1 #K
 
31
#alpha = 1 #/K
 
32
 
 
33
def wavenumber():
 
34
  wavelength = wavelengthfactor*D
 
35
  return 2.0*pi/wavelength
 
36
 
 
37
def nond_wavenumber():
 
38
  wavelength = wavelengthfactor
 
39
  return 2.0*pi/wavelength
 
40
 
 
41
def t0_eta():
 
42
  delta_rho = (rho0-rhou)*g
 
43
  rhog = rho0*g
 
44
  mu = mu0
 
45
  k = wavenumber()
 
46
  return -1.0*(((delta_rho*k**2*mu - k**2*mu*rhog)*D*sinh(D*k)**2 - (delta_rho*k**2*mu - k**2*mu*rhog)*D*cosh(D*k)**2 - (delta_rho*k*mu -\
 
47
k*mu*rhog)*sinh(D*k)*cosh(D*k) - sqrt((delta_rho**2*k**2 - 2*delta_rho*k**2*rhog + k**2*rhog**2)*D**2*cosh(D*k)**4 - 2*(delta_rho**2*k -\
 
48
2*delta_rho*k*rhog + k*rhog**2)*D*sinh(D*k)**3*cosh(D*k) + 2*(delta_rho**2*k - 2*delta_rho*k*rhog + k*rhog**2)*D*sinh(D*k)*cosh(D*k)**3 +\
 
49
((delta_rho**2*k**2 + 2*delta_rho*k**2*rhog + k**2*rhog**2)*D**2 + 4*delta_rho*rhog)*sinh(D*k)**4 - (2*(delta_rho**2*k**2 + k**2*rhog**2)*D**2 -\
 
50
delta_rho**2 + 2*delta_rho*rhog - rhog**2)*sinh(D*k)**2*cosh(D*k)**2)*k*mu)/(delta_rho*rhog*sinh(D*k)**2))
 
51
 
 
52
def t0_xi():
 
53
  delta_rho = (rho0-rhou)*g
 
54
  rhog = rho0*g
 
55
  mu = mu0
 
56
  k = wavenumber()
 
57
  return -1.0*(((delta_rho*k**2*mu - k**2*mu*rhog)*D*sinh(D*k)**2 - (delta_rho*k**2*mu - k**2*mu*rhog)*D*cosh(D*k)**2 - (delta_rho*k*mu -\
 
58
k*mu*rhog)*sinh(D*k)*cosh(D*k) + sqrt((delta_rho**2*k**2 - 2*delta_rho*k**2*rhog + k**2*rhog**2)*D**2*cosh(D*k)**4 - 2*(delta_rho**2*k -\
 
59
2*delta_rho*k*rhog + k*rhog**2)*D*sinh(D*k)**3*cosh(D*k) + 2*(delta_rho**2*k - 2*delta_rho*k*rhog + k*rhog**2)*D*sinh(D*k)*cosh(D*k)**3 +\
 
60
((delta_rho**2*k**2 + 2*delta_rho*k**2*rhog + k**2*rhog**2)*D**2 + 4*delta_rho*rhog)*sinh(D*k)**4 - (2*(delta_rho**2*k**2 + k**2*rhog**2)*D**2 -\
 
61
delta_rho**2 + 2*delta_rho*rhog - rhog**2)*sinh(D*k)**2*cosh(D*k)**2)*k*mu)/(delta_rho*rhog*sinh(D*k)**2))
 
62
 
 
63
def t0():
 
64
  return min(t0_eta(), t0_xi())
 
65
 
 
66
def nond_factor():
 
67
  return rho0*g*D*t0()/mu0
 
68
 
 
69
def nond_eta0():
 
70
  return eta0/D
 
71
 
 
72
def nond_xi0():
 
73
  return xi0/D
 
74
 
 
75
def nond_F(x,t):
 
76
  k = nond_wavenumber()
 
77
  F0 = nond_eta0()
 
78
  G0 = nond_xi0()
 
79
  delta_rho = (rho0-rhou)*nond_factor()/rho0
 
80
  zprime = depth/D
 
81
  T0 = deltaT
 
82
  alphag = nond_factor()*alpha*T0
 
83
  rhog = nond_factor() # use this as a proxy for the nondimensional factorisation
 
84
  return ((-0.5*(exp(-delta_rho*rhog*t*sinh(k)**2/((delta_rho*k - k*rhog)*sinh(k)*cosh(k) + delta_rho*k**2 - k**2*rhog -\
 
85
sqrt((delta_rho**2 + 2*delta_rho*rhog + rhog**2)*sinh(k)**4 + delta_rho**2*k**2 - 2*delta_rho*k**2*rhog + k**2*rhog**2 + 2*(delta_rho**2*k -\
 
86
2*delta_rho*k*rhog + k*rhog**2)*sinh(k)*cosh(k) - (4*delta_rho*k**2*rhog - delta_rho**2 + 2*delta_rho*rhog -\
 
87
rhog**2)*sinh(k)**2)*k))/((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 - rhog*sinh(k)*cosh(k)**2)/((delta_rho +\
 
88
rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2 - sqrt(delta_rho**2*k**2*sinh(k)**4 -\
 
89
2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4\
 
90
+ k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) +\
 
91
2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k) - 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 -\
 
92
2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 -\
 
93
2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 + rhog**2*sinh(k)**2*cosh(k)**2)) - (k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3\
 
94
- rhog*sinh(k)*cosh(k)**2)/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2\
 
95
  + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
96
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
97
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
98
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
99
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
100
rhog**2*sinh(k)**2*cosh(k)**2))) - exp(-delta_rho*rhog*t*sinh(k)**2/((delta_rho*k - k*rhog)*sinh(k)*cosh(k) + delta_rho*k**2 - k**2*rhog +\
 
101
sqrt((delta_rho**2 + 2*delta_rho*rhog + rhog**2)*sinh(k)**4 + delta_rho**2*k**2 - 2*delta_rho*k**2*rhog + k**2*rhog**2 + 2*(delta_rho**2*k -\
 
102
2*delta_rho*k*rhog + k*rhog**2)*sinh(k)*cosh(k) - (4*delta_rho*k**2*rhog - delta_rho**2 + 2*delta_rho*rhog -\
 
103
rhog**2)*sinh(k)**2)*k))/((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 - rhog*sinh(k)*cosh(k)**2)/((delta_rho +\
 
104
rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2 - sqrt(delta_rho**2*k**2*sinh(k)**4 -\
 
105
2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4\
 
106
+ k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) +\
 
107
2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k) - 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 -\
 
108
2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 -\
 
109
2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 + rhog**2*sinh(k)**2*cosh(k)**2)) - (k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3\
 
110
- rhog*sinh(k)*cosh(k)**2)/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2\
 
111
  + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
112
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
113
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
114
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
115
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
116
rhog**2*sinh(k)**2*cosh(k)**2))))*(G0 + (alphag*k*sinh(k*zprime)*cosh(k) - (alphag*k*zprime*cosh(k*zprime) -\
 
117
alphag*sinh(k*zprime))*sinh(k))/(T0*delta_rho*sinh(k)**2)) + (F0 - ((alphag*k*zprime*cosh(k*zprime) -\
 
118
alphag*sinh(k*zprime))*sinh(k)*cosh(k) - (alphag*k*zprime*sinh(k*zprime) - alphag*cosh(k*zprime))*sinh(k)**2 -\
 
119
alphag*k*sinh(k*zprime))/(T0*rhog*sinh(k)**2))*(((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 -\
 
120
rhog*sinh(k)*cosh(k)**2)/(((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 - rhog*sinh(k)*cosh(k)**2)/((delta_rho +\
 
121
rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2 - sqrt(delta_rho**2*k**2*sinh(k)**4 -\
 
122
2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4\
 
123
+ k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) +\
 
124
2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k) - 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 -\
 
125
2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 -\
 
126
2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 + rhog**2*sinh(k)**2*cosh(k)**2)) - (k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3\
 
127
- rhog*sinh(k)*cosh(k)**2)/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2\
 
128
  + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
129
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
130
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
131
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
132
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
133
rhog**2*sinh(k)**2*cosh(k)**2)))*((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k +\
 
134
k*rhog)*cosh(k)**2 + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
135
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
136
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
137
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
138
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
139
rhog**2*sinh(k)**2*cosh(k)**2))) + 1)*exp(-delta_rho*rhog*t*sinh(k)**2/((delta_rho*k - k*rhog)*sinh(k)*cosh(k) + delta_rho*k**2 - k**2*rhog\
 
140
+ sqrt((delta_rho**2 + 2*delta_rho*rhog + rhog**2)*sinh(k)**4 + delta_rho**2*k**2 - 2*delta_rho*k**2*rhog + k**2*rhog**2 + 2*(delta_rho**2*k\
 
141
- 2*delta_rho*k*rhog + k*rhog**2)*sinh(k)*cosh(k) - (4*delta_rho*k**2*rhog - delta_rho**2 + 2*delta_rho*rhog - rhog**2)*sinh(k)**2)*k)) -\
 
142
  (k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 -\
 
143
rhog*sinh(k)*cosh(k)**2)*exp(-delta_rho*rhog*t*sinh(k)**2/((delta_rho*k - k*rhog)*sinh(k)*cosh(k) + delta_rho*k**2 - k**2*rhog -\
 
144
sqrt((delta_rho**2 + 2*delta_rho*rhog + rhog**2)*sinh(k)**4 + delta_rho**2*k**2 - 2*delta_rho*k**2*rhog + k**2*rhog**2 + 2*(delta_rho**2*k -\
 
145
2*delta_rho*k*rhog + k*rhog**2)*sinh(k)*cosh(k) - (4*delta_rho*k**2*rhog - delta_rho**2 + 2*delta_rho*rhog -\
 
146
rhog**2)*sinh(k)**2)*k))/(((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 - rhog*sinh(k)*cosh(k)**2)/((delta_rho +\
 
147
rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2 - sqrt(delta_rho**2*k**2*sinh(k)**4 -\
 
148
2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4\
 
149
+ k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) +\
 
150
2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k) - 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 -\
 
151
2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 -\
 
152
2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 + rhog**2*sinh(k)**2*cosh(k)**2)) - (k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3\
 
153
- rhog*sinh(k)*cosh(k)**2)/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2\
 
154
  + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
155
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
156
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
157
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
158
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
159
rhog**2*sinh(k)**2*cosh(k)**2)))*((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k +\
 
160
k*rhog)*cosh(k)**2 + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
161
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
162
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
163
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
164
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
165
rhog**2*sinh(k)**2*cosh(k)**2)))) + ((alphag*k*zprime*cosh(k*zprime) - alphag*sinh(k*zprime))*sinh(k)*cosh(k) -\
 
166
(alphag*k*zprime*sinh(k*zprime) - alphag*cosh(k*zprime))*sinh(k)**2 - alphag*k*sinh(k*zprime))/(T0*rhog*sinh(k)**2)))*cos(k*x)
 
167
 
 
168
def nond_G(x,t):
 
169
  k = nond_wavenumber()
 
170
  F0 = nond_eta0()
 
171
  G0 = nond_xi0()
 
172
  delta_rho = (rho0-rhou)*nond_factor()/rho0
 
173
  zprime = depth/D
 
174
  T0 = deltaT
 
175
  alphag = nond_factor()*alpha*T0
 
176
  rhog = nond_factor() # use this as a proxy for the nondimensional factorisation
 
177
  return ((((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 -\
 
178
rhog*sinh(k)*cosh(k)**2)*exp(-delta_rho*rhog*t*sinh(k)**2/((delta_rho*k - k*rhog)*sinh(k)*cosh(k) + delta_rho*k**2 - k**2*rhog -\
 
179
sqrt((delta_rho**2 + 2*delta_rho*rhog + rhog**2)*sinh(k)**4 + delta_rho**2*k**2 - 2*delta_rho*k**2*rhog + k**2*rhog**2 + 2*(delta_rho**2*k -\
 
180
2*delta_rho*k*rhog + k*rhog**2)*sinh(k)*cosh(k) - (4*delta_rho*k**2*rhog - delta_rho**2 + 2*delta_rho*rhog -\
 
181
rhog**2)*sinh(k)**2)*k))/(((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 - rhog*sinh(k)*cosh(k)**2)/((delta_rho +\
 
182
rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2 - sqrt(delta_rho**2*k**2*sinh(k)**4 -\
 
183
2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4\
 
184
+ k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) +\
 
185
2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k) - 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 -\
 
186
2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 -\
 
187
2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 + rhog**2*sinh(k)**2*cosh(k)**2)) - (k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3\
 
188
- rhog*sinh(k)*cosh(k)**2)/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2\
 
189
  + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
190
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
191
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
192
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
193
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
194
rhog**2*sinh(k)**2*cosh(k)**2)))*((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k +\
 
195
k*rhog)*cosh(k)**2 - sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
196
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
197
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
198
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
199
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
200
rhog**2*sinh(k)**2*cosh(k)**2))) - (k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 -\
 
201
rhog*sinh(k)*cosh(k)**2)*exp(-delta_rho*rhog*t*sinh(k)**2/((delta_rho*k - k*rhog)*sinh(k)*cosh(k) + delta_rho*k**2 - k**2*rhog +\
 
202
sqrt((delta_rho**2 + 2*delta_rho*rhog + rhog**2)*sinh(k)**4 + delta_rho**2*k**2 - 2*delta_rho*k**2*rhog + k**2*rhog**2 + 2*(delta_rho**2*k -\
 
203
2*delta_rho*k*rhog + k*rhog**2)*sinh(k)*cosh(k) - (4*delta_rho*k**2*rhog - delta_rho**2 + 2*delta_rho*rhog -\
 
204
rhog**2)*sinh(k)**2)*k))/(((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 - rhog*sinh(k)*cosh(k)**2)/((delta_rho +\
 
205
rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2 - sqrt(delta_rho**2*k**2*sinh(k)**4 -\
 
206
2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4\
 
207
+ k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) +\
 
208
2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k) - 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 -\
 
209
2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 -\
 
210
2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 + rhog**2*sinh(k)**2*cosh(k)**2)) - (k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3\
 
211
- rhog*sinh(k)*cosh(k)**2)/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2\
 
212
  + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
213
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
214
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
215
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
216
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
217
rhog**2*sinh(k)**2*cosh(k)**2)))*((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k +\
 
218
k*rhog)*cosh(k)**2 + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
219
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
220
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
221
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
222
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
223
rhog**2*sinh(k)**2*cosh(k)**2))))*(G0 + (alphag*k*sinh(k*zprime)*cosh(k) - (alphag*k*zprime*cosh(k*zprime) -\
 
224
alphag*sinh(k*zprime))*sinh(k))/(T0*delta_rho*sinh(k)**2)) - 2*(F0 - ((alphag*k*zprime*cosh(k*zprime) -\
 
225
alphag*sinh(k*zprime))*sinh(k)*cosh(k) - (alphag*k*zprime*sinh(k*zprime) - alphag*cosh(k*zprime))*sinh(k)**2 -\
 
226
alphag*k*sinh(k*zprime))/(T0*rhog*sinh(k)**2))*(((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 -\
 
227
rhog*sinh(k)*cosh(k)**2)/(((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 - rhog*sinh(k)*cosh(k)**2)/((delta_rho +\
 
228
rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2 - sqrt(delta_rho**2*k**2*sinh(k)**4 -\
 
229
2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4\
 
230
+ k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) +\
 
231
2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k) - 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 -\
 
232
2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 -\
 
233
2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 + rhog**2*sinh(k)**2*cosh(k)**2)) - (k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3\
 
234
- rhog*sinh(k)*cosh(k)**2)/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2\
 
235
  + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
236
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
237
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
238
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
239
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
240
rhog**2*sinh(k)**2*cosh(k)**2)))*((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k +\
 
241
k*rhog)*cosh(k)**2 + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
242
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
243
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
244
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
245
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
246
rhog**2*sinh(k)**2*cosh(k)**2))) + 1)*(k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 -\
 
247
rhog*sinh(k)*cosh(k)**2)*exp(-delta_rho*rhog*t*sinh(k)**2/((delta_rho*k - k*rhog)*sinh(k)*cosh(k) + delta_rho*k**2 - k**2*rhog +\
 
248
sqrt((delta_rho**2 + 2*delta_rho*rhog + rhog**2)*sinh(k)**4 + delta_rho**2*k**2 - 2*delta_rho*k**2*rhog + k**2*rhog**2 + 2*(delta_rho**2*k -\
 
249
2*delta_rho*k*rhog + k*rhog**2)*sinh(k)*cosh(k) - (4*delta_rho*k**2*rhog - delta_rho**2 + 2*delta_rho*rhog -\
 
250
rhog**2)*sinh(k)**2)*k))/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2 +\
 
251
sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4\
 
252
- 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 -\
 
253
  2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k) -\
 
254
4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2\
 
255
+ 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 + rhog**2*sinh(k)**2*cosh(k)**2)) - (k*rhog*sinh(k)**2*cosh(k) -\
 
256
k*rhog*cosh(k)**3 + rhog*sinh(k)**3 - rhog*sinh(k)*cosh(k)**2)**2*exp(-delta_rho*rhog*t*sinh(k)**2/((delta_rho*k - k*rhog)*sinh(k)*cosh(k)\
 
257
+ delta_rho*k**2 - k**2*rhog - sqrt((delta_rho**2 + 2*delta_rho*rhog + rhog**2)*sinh(k)**4 + delta_rho**2*k**2 - 2*delta_rho*k**2*rhog +\
 
258
k**2*rhog**2 + 2*(delta_rho**2*k - 2*delta_rho*k*rhog + k*rhog**2)*sinh(k)*cosh(k) - (4*delta_rho*k**2*rhog - delta_rho**2 +\
 
259
2*delta_rho*rhog - rhog**2)*sinh(k)**2)*k))/(((k*rhog*sinh(k)**2*cosh(k) - k*rhog*cosh(k)**3 + rhog*sinh(k)**3 -\
 
260
rhog*sinh(k)*cosh(k)**2)/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k + k*rhog)*cosh(k)**2 -\
 
261
sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4\
 
262
- 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 -\
 
263
  2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k) -\
 
264
4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2\
 
265
+ 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 + rhog**2*sinh(k)**2*cosh(k)**2)) - (k*rhog*sinh(k)**2*cosh(k) -\
 
266
k*rhog*cosh(k)**3 + rhog*sinh(k)**3 - rhog*sinh(k)*cosh(k)**2)/((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 +\
 
267
(delta_rho*k + k*rhog)*cosh(k)**2 + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 +\
 
268
delta_rho**2*k**2*cosh(k)**4 + 2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 -\
 
269
2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 + k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 +\
 
270
4*delta_rho*k*rhog*sinh(k)**3*cosh(k) - 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) +\
 
271
2*k*rhog**2*sinh(k)*cosh(k)**3 + delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
272
rhog**2*sinh(k)**2*cosh(k)**2)))*((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k +\
 
273
k*rhog)*cosh(k)**2 - sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
274
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
275
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
276
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
277
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
278
rhog**2*sinh(k)**2*cosh(k)**2))*((delta_rho + rhog)*sinh(k)*cosh(k) - (delta_rho*k + k*rhog)*sinh(k)**2 + (delta_rho*k +\
 
279
k*rhog)*cosh(k)**2 + sqrt(delta_rho**2*k**2*sinh(k)**4 - 2*delta_rho**2*k**2*sinh(k)**2*cosh(k)**2 + delta_rho**2*k**2*cosh(k)**4 +\
 
280
2*delta_rho*k**2*rhog*sinh(k)**4 - 2*delta_rho*k**2*rhog*cosh(k)**4 + k**2*rhog**2*sinh(k)**4 - 2*k**2*rhog**2*sinh(k)**2*cosh(k)**2 +\
 
281
k**2*rhog**2*cosh(k)**4 - 2*delta_rho**2*k*sinh(k)**3*cosh(k) + 2*delta_rho**2*k*sinh(k)*cosh(k)**3 + 4*delta_rho*k*rhog*sinh(k)**3*cosh(k)\
 
282
- 4*delta_rho*k*rhog*sinh(k)*cosh(k)**3 - 2*k*rhog**2*sinh(k)**3*cosh(k) + 2*k*rhog**2*sinh(k)*cosh(k)**3 +\
 
283
  delta_rho**2*sinh(k)**2*cosh(k)**2 + 4*delta_rho*rhog*sinh(k)**4 - 2*delta_rho*rhog*sinh(k)**2*cosh(k)**2 +\
 
284
rhog**2*sinh(k)**2*cosh(k)**2)))) - (alphag*k*sinh(k*zprime)*cosh(k) - (alphag*k*zprime*cosh(k*zprime) -\
 
285
alphag*sinh(k*zprime))*sinh(k))/(T0*delta_rho*sinh(k)**2)))*cos(k*x)
 
286
 
 
287
def nond_amp(t):
 
288
  return max(nond_F(0.0,t), nond_G(0.0,t))
 
289
 
 
290
def numerical_amp(statn, t):
 
291
  return interp([t], statn["ElapsedTime"]["value"], statn["Fluid"]["FreeSurface"]["max"])[0]
 
292
 
 
293
def nond_error_amp(statn, t):
 
294
  return abs(nond_amp(t)-numerical_amp(statn, t))
 
295
 
 
296
def nond_F_ss(x):
 
297
  k = nond_wavenumber()
 
298
  F0 = nond_eta0()
 
299
  G0 = nond_xi0()
 
300
  delta_rho = (rho0-rhou)*nond_factor()/rho0
 
301
  zprime = depth/D
 
302
  T0 = deltaT
 
303
  alphag = nond_factor()
 
304
  rhog = nond_factor() # use this as a proxy for the nondimensional factorisation
 
305
  return (((alphag*k*zprime*cosh(k*zprime) - alphag*sinh(k*zprime))*sinh(k)*cosh(k) - (alphag*k*zprime*sinh(k*zprime) - \
 
306
alphag*cosh(k*zprime))*sinh(k)**2 - alphag*k*sinh(k*zprime))/(T0*rhog*sinh(k)**2))*cos(k*x)
 
307
 
 
308
 
 
309
def nond_G_ss(x):
 
310
  k = nond_wavenumber()
 
311
  F0 = nond_eta0()
 
312
  G0 = nond_xi0()
 
313
  delta_rho = (rho0-rhou)*nond_factor()/rho0
 
314
  zprime = depth/D
 
315
  T0 = deltaT
 
316
  alphag = nond_factor()
 
317
  rhog = nond_factor() # use this as a proxy for the nondimensional factorisation
 
318
  return (-(alphag*k*sinh(k*zprime)*cosh(k) - (alphag*k*zprime*cosh(k*zprime) - alphag*sinh(k*zprime))*sinh(k))/\
 
319
(T0*delta_rho*sinh(k)**2))*cos(k*x)
 
320
 
 
321