A Modular Arbitrary-Order Ocean-Atmosphere Model -- Stochastic implementation
Functions/Subroutines | Variables
corrmod Module Reference

Module to initialize the correlation matrix of the unresolved variables. More...

Functions/Subroutines

subroutine, public init_corr
 Subroutine to initialise the computation of the correlation. More...
 
subroutine corrcomp_from_def (s)
 Subroutine to compute the correlation of the unresolved variables \(\langle Y \otimes Y^s \rangle\) at time \(s\) from the definition given inside the module. More...
 
subroutine corrcomp_from_spline (s)
 Subroutine to compute the correlation of the unresolved variables \(\langle Y \otimes Y^s \rangle\) at time \(s\) from the spline representation. More...
 
subroutine splint (xa, ya, y2a, n, x, y)
 Routine to compute the spline representation parameters. More...
 
real(kind=8) function fs (s, p)
 Exponential fit function. More...
 
subroutine corrcomp_from_fit (s)
 Subroutine to compute the correlation of the unresolved variables \(\langle Y \otimes Y^s \rangle\) at time \(s\) from the exponential representation. More...
 

Variables

real(kind=8), dimension(:), allocatable, public mean
 Vector holding the mean of the unresolved dynamics (reduced version) More...
 
real(kind=8), dimension(:), allocatable, public mean_full
 Vector holding the mean of the unresolved dynamics (full version) More...
 
real(kind=8), dimension(:,:), allocatable, public corr_i_full
 Covariance matrix of the unresolved variables (full version) More...
 
real(kind=8), dimension(:,:), allocatable, public inv_corr_i_full
 Inverse of the covariance matrix of the unresolved variables (full version) More...
 
real(kind=8), dimension(:,:), allocatable, public corr_i
 Covariance matrix of the unresolved variables (reduced version) More...
 
real(kind=8), dimension(:,:), allocatable, public inv_corr_i
 Inverse of the covariance matrix of the unresolved variables (reduced version) More...
 
real(kind=8), dimension(:,:), allocatable, public corr_ij
 Matrix holding the correlation matrix at a given time. More...
 
real(kind=8), dimension(:,:,:), allocatable y2
 Vector holding coefficient of the spline and exponential correlation representation. More...
 
real(kind=8), dimension(:,:,:), allocatable ya
 Vector holding coefficient of the spline and exponential correlation representation. More...
 
real(kind=8), dimension(:), allocatable xa
 Vector holding coefficient of the spline and exponential correlation representation. More...
 
integer nspl
 Integers needed by the spline representation of the correlation. More...
 
integer klo
 
integer khi
 
procedure(corrcomp_from_spline), pointer, public corrcomp
 Pointer to the correlation computation routine. More...
 

Detailed Description

Module to initialize the correlation matrix of the unresolved variables.

Remarks

Function/Subroutine Documentation

subroutine corrmod::corrcomp_from_def ( real(kind=8), intent(in)  s)
private

Subroutine to compute the correlation of the unresolved variables \(\langle Y \otimes Y^s \rangle\) at time \(s\) from the definition given inside the module.

Parameters
stime \(s\) at which the correlation is computed

Definition at line 148 of file corrmod.f90.

148  REAL(KIND=8), INTENT(IN) :: s
149  REAL(KIND=8) :: y
150  INTEGER :: i,j
151 
152  ! Definition of the corr_ij matrix as a function of time
153 
154  corr_ij(10,10)=((7.66977252307669*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (1.0240906173830213*cos(&
155  &0.07283568782600224*s))/exp(0.017262015588746404*s) - (0.6434985372062336*sin(0.039597160512071454*s&
156  &))/exp(0.06567483898489704*s) + (0.6434985372062335*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
157  corr_ij(10,9)=((0.6434985372062321*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) - (0.6434985372062324*co&
158  &s(0.07283568782600224*s))/exp(0.017262015588746404*s) + (7.669772523076694*sin(0.039597160512071454*&
159  &s))/exp(0.06567483898489704*s) + (1.024090617383021*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
160  corr_ij(10,8)=0
161  corr_ij(10,7)=0
162  corr_ij(10,6)=0
163  corr_ij(10,5)=((-2.236364132659011*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (6.952804148086198*cos&
164  &(0.07283568782600224*s))/exp(0.017262015588746404*s) - (1.4494534432272481*sin(0.039597160512071454*&
165  &s))/exp(0.06567483898489704*s) - (0.6818177416446283*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
166  corr_ij(10,4)=((1.4494534432272483*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (0.6818177416446293*co&
167  &s(0.07283568782600224*s))/exp(0.017262015588746404*s) - (2.2363641326590127*sin(0.039597160512071454&
168  &*s))/exp(0.06567483898489704*s) + (6.952804148086195*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
169  corr_ij(10,3)=0
170  corr_ij(10,2)=0
171  corr_ij(10,1)=0
172  corr_ij(9,10)=((-0.6434985372062344*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (0.643498537206234*co&
173  &s(0.07283568782600224*s))/exp(0.017262015588746404*s) - (7.669772523076689*sin(0.039597160512071454*&
174  &s))/exp(0.06567483898489704*s) - (1.0240906173830204*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
175  corr_ij(9,9)=((7.66977252307669*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (1.0240906173830204*cos(&
176  &0.07283568782600224*s))/exp(0.017262015588746404*s) - (0.643498537206233*sin(0.039597160512071454*s)&
177  &)/exp(0.06567483898489704*s) + (0.6434985372062327*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
178  corr_ij(9,8)=0
179  corr_ij(9,7)=0
180  corr_ij(9,6)=0
181  corr_ij(9,5)=((-1.4494534432272477*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) - (0.6818177416446249*c&
182  &os(0.07283568782600224*s))/exp(0.017262015588746404*s) + (2.2363641326590105*sin(0.03959716051207145&
183  &4*s))/exp(0.06567483898489704*s) - (6.952804148086195*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
184  corr_ij(9,4)=((-2.2363641326590127*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (6.952804148086194*co&
185  &s(0.07283568782600224*s))/exp(0.017262015588746404*s) - (1.4494534432272486*sin(0.039597160512071454&
186  &*s))/exp(0.06567483898489704*s) - (0.6818177416446249*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
187  corr_ij(9,3)=0
188  corr_ij(9,2)=0
189  corr_ij(9,1)=0
190  corr_ij(8,10)=0
191 
192  corr_ij(8,9)=0
193  corr_ij(8,8)=(9.135647293470983/exp(0.05076718239027029*s) + 2.2233889637758932/exp(0.016285467000648854*s))
194  corr_ij(8,7)=0
195  corr_ij(8,6)=0
196  corr_ij(8,5)=0
197  corr_ij(8,4)=0
198  corr_ij(8,3)=(-5.938566084855411/exp(0.05076718239027029*s) + 11.97129741027622/exp(0.016285467000648854*s))
199  corr_ij(8,2)=0
200  corr_ij(8,1)=0
201  corr_ij(7,10)=0
202 
203  corr_ij(7,9)=0
204  corr_ij(7,8)=0
205  corr_ij(7,7)=((11.518026982819887*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (0.05351107793747755*c&
206  &os(0.11425932545092894*s))/exp(0.019700737327669783*s) - (0.14054811601869432*sin(0.0293414097268719&
207  &26*s))/exp(0.04435489221745234*s) + (0.14054811601869702*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
208  corr_ij(7,6)=((0.14054811601869532*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) - (0.14054811601869702*&
209  &cos(0.11425932545092894*s))/exp(0.019700737327669783*s) + (11.518026982819887*sin(0.0293414097268719&
210  &26*s))/exp(0.04435489221745234*s) + (0.0535110779374777*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
211  corr_ij(7,5)=0
212  corr_ij(7,4)=0
213  corr_ij(7,3)=0
214  corr_ij(7,2)=((-0.732907009016115*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (2.728845031386875*cos&
215  &(0.11425932545092894*s))/exp(0.019700737327669783*s) - (2.4717920234033532*sin(0.029341409726871926*&
216  &s))/exp(0.04435489221745234*s) - (0.24003801347124257*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
217  corr_ij(7,1)=((2.4717920234033532*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (0.2400380134712426*co&
218  &s(0.11425932545092894*s))/exp(0.019700737327669783*s) - (0.7329070090161153*sin(0.029341409726871926&
219  &*s))/exp(0.04435489221745234*s) + (2.728845031386876*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
220  corr_ij(6,10)=0
221 
222  corr_ij(6,9)=0
223  corr_ij(6,8)=0
224  corr_ij(6,7)=((-0.1405481160186977*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (0.14054811601869713*&
225  &cos(0.11425932545092894*s))/exp(0.019700737327669783*s) - (11.518026982819885*sin(0.0293414097268719&
226  &26*s))/exp(0.04435489221745234*s) - (0.05351107793747755*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
227  corr_ij(6,6)=((11.518026982819885*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (0.05351107793747768*c&
228  &os(0.11425932545092894*s))/exp(0.019700737327669783*s) - (0.14054811601869832*sin(0.0293414097268719&
229  &26*s))/exp(0.04435489221745234*s) + (0.14054811601869707*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
230  corr_ij(6,5)=0
231  corr_ij(6,4)=0
232  corr_ij(6,3)=0
233  corr_ij(6,2)=((-2.471792023403353*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) - (0.2400380134712425*co&
234  &s(0.11425932545092894*s))/exp(0.019700737327669783*s) + (0.7329070090161155*sin(0.029341409726871926&
235  &*s))/exp(0.04435489221745234*s) - (2.7288450313868755*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
236  corr_ij(6,1)=((-0.7329070090161154*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (2.728845031386876*co&
237  &s(0.11425932545092894*s))/exp(0.019700737327669783*s) - (2.4717920234033524*sin(0.029341409726871926&
238  &*s))/exp(0.04435489221745234*s) - (0.24003801347124343*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
239  corr_ij(5,10)=((0.5794534449999711*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (4.136986570427212*cos&
240  &(0.07283568782600224*s))/exp(0.017262015588746404*s) - (1.0360597341248128*sin(0.039597160512071454*&
241  &s))/exp(0.06567483898489704*s) + (3.167330918996692*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
242  corr_ij(5,9)=((1.0360597341248134*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) - (3.1673309189966856*co&
243  &s(0.07283568782600224*s))/exp(0.017262015588746404*s) + (0.5794534449999746*sin(0.039597160512071454&
244  &*s))/exp(0.06567483898489704*s) + (4.1369865704272115*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
245  corr_ij(5,8)=0
246  corr_ij(5,7)=0
247  corr_ij(5,6)=0
248  corr_ij(5,5)=((-0.37825091063447547*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (30.094690926061638*&
249  &cos(0.07283568782600224*s))/exp(0.017262015588746404*s) + (0.16085380971100194*sin(0.039597160512071&
250  &454*s))/exp(0.06567483898489704*s) - (0.1608538097109995*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
251  corr_ij(5,4)=((-0.16085380971100238*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (0.16085380971100127&
252  &*cos(0.07283568782600224*s))/exp(0.017262015588746404*s) - (0.37825091063447586*sin(0.03959716051207&
253  &1454*s))/exp(0.06567483898489704*s) + (30.09469092606163*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
254  corr_ij(5,3)=0
255  corr_ij(5,2)=0
256  corr_ij(5,1)=0
257  corr_ij(4,10)=((-1.0360597341248106*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (3.167330918996689*co&
258  &s(0.07283568782600224*s))/exp(0.017262015588746404*s) - (0.5794534449999716*sin(0.039597160512071454&
259  &*s))/exp(0.06567483898489704*s) - (4.1369865704272115*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
260  corr_ij(4,9)=((0.5794534449999711*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (4.1369865704272115*co&
261  &s(0.07283568782600224*s))/exp(0.017262015588746404*s) - (1.0360597341248114*sin(0.039597160512071454&
262  &*s))/exp(0.06567483898489704*s) + (3.1673309189966843*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
263  corr_ij(4,8)=0
264  corr_ij(4,7)=0
265  corr_ij(4,6)=0
266  corr_ij(4,5)=((0.16085380971100194*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) - (0.16085380971100371*&
267  &cos(0.07283568782600224*s))/exp(0.017262015588746404*s) + (0.37825091063447497*sin(0.039597160512071&
268  &454*s))/exp(0.06567483898489704*s) - (30.094690926061617*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
269  corr_ij(4,4)=((-0.37825091063447536*cos(0.039597160512071454*s))/exp(0.06567483898489704*s) + (30.094690926061617*&
270  &cos(0.07283568782600224*s))/exp(0.017262015588746404*s) + (0.16085380971100172*sin(0.039597160512071&
271  &454*s))/exp(0.06567483898489704*s) - (0.16085380971100616*sin(0.07283568782600224*s))/exp(0.017262015588746404*s))
272  corr_ij(4,3)=0
273  corr_ij(4,2)=0
274  corr_ij(4,1)=0
275  corr_ij(3,10)=0
276 
277  corr_ij(3,9)=0
278  corr_ij(3,8)=(0.24013456462471527/exp(0.05076718239027029*s) + 5.792596760796093/exp(0.016285467000648854*s))
279  corr_ij(3,7)=0
280  corr_ij(3,6)=0
281  corr_ij(3,5)=0
282  corr_ij(3,4)=0
283  corr_ij(3,3)=(-0.15609785880208227/exp(0.05076718239027029*s) + 31.18882918422289/exp(0.016285467000648854*s))
284  corr_ij(3,2)=0
285  corr_ij(3,1)=0
286  corr_ij(2,10)=0
287 
288  corr_ij(2,9)=0
289  corr_ij(2,8)=0
290  corr_ij(2,7)=((1.6172201305728584*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (0.37871789179790255*c&
291  &os(0.11425932545092894*s))/exp(0.019700737327669783*s) + (1.2889451151208258*sin(0.02934140972687192&
292  &6*s))/exp(0.04435489221745234*s) + (1.4228849217537705*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
293  corr_ij(2,6)=((-1.2889451151208255*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) - (1.4228849217537702*c&
294  &os(0.11425932545092894*s))/exp(0.019700737327669783*s) + (1.6172201305728586*sin(0.02934140972687192&
295  &6*s))/exp(0.04435489221745234*s) + (0.3787178917979035*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
296  corr_ij(2,5)=0
297  corr_ij(2,4)=0
298  corr_ij(2,3)=0
299  corr_ij(2,2)=((0.1789135645266575*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (26.817024457844113*co&
300  &s(0.11425932545092894*s))/exp(0.019700737327669783*s) - (0.4268927977731004*sin(0.029341409726871926&
301  &*s))/exp(0.04435489221745234*s) + (0.4268927977730982*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
302  corr_ij(2,1)=((0.4268927977731007*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) - (0.42689279777309963*c&
303  &os(0.11425932545092894*s))/exp(0.019700737327669783*s) + (0.17891356452665746*sin(0.0293414097268719&
304  &26*s))/exp(0.04435489221745234*s) + (26.81702445784412*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
305  corr_ij(1,10)=0
306 
307  corr_ij(1,9)=0
308  corr_ij(1,8)=0
309  corr_ij(1,7)=((1.288945115120824*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (1.4228849217537711*cos&
310  &(0.11425932545092894*s))/exp(0.019700737327669783*s) - (1.617220130572856*sin(0.029341409726871926*s&
311  &))/exp(0.04435489221745234*s) - (0.3787178917979028*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
312  corr_ij(1,6)=((1.6172201305728564*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (0.37871789179790377*c&
313  &os(0.11425932545092894*s))/exp(0.019700737327669783*s) + (1.2889451151208242*sin(0.02934140972687192&
314  &6*s))/exp(0.04435489221745234*s) + (1.4228849217537711*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
315  corr_ij(1,5)=0
316  corr_ij(1,4)=0
317  corr_ij(1,3)=0
318  corr_ij(1,2)=((-0.4268927977731002*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (0.4268927977730981*c&
319  &os(0.11425932545092894*s))/exp(0.019700737327669783*s) - (0.1789135645266573*sin(0.02934140972687192&
320  &6*s))/exp(0.04435489221745234*s) - (26.81702445784412*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
321  corr_ij(1,1)=((0.1789135645266574*cos(0.029341409726871926*s))/exp(0.04435489221745234*s) + (26.817024457844113*co&
322  &s(0.11425932545092894*s))/exp(0.019700737327669783*s) - (0.42689279777310024*sin(0.02934140972687192&
323  &6*s))/exp(0.04435489221745234*s) + (0.4268927977730997*sin(0.11425932545092894*s))/exp(0.019700737327669783*s))
324 
325  corr_ij=q_au**2*corr_ij
326 
subroutine corrmod::corrcomp_from_fit ( real(kind=8), intent(in)  s)
private

Subroutine to compute the correlation of the unresolved variables \(\langle Y \otimes Y^s \rangle\) at time \(s\) from the exponential representation.

Parameters
stime \(s\) at which the correlation is computed

Definition at line 399 of file corrmod.f90.

399  REAL(KIND=8), INTENT(IN) :: s
400  REAL(KIND=8) :: y
401  INTEGER :: i,j
402 
403  corr_ij=0.d0
404  DO i=1,n_unres
405  DO j=1,n_unres
406  corr_ij(i,j)=fs(s,ya(i,j,:))
407  END DO
408  END DO
subroutine corrmod::corrcomp_from_spline ( real(kind=8), intent(in)  s)
private

Subroutine to compute the correlation of the unresolved variables \(\langle Y \otimes Y^s \rangle\) at time \(s\) from the spline representation.

Parameters
stime \(s\) at which the correlation is computed

Definition at line 333 of file corrmod.f90.

333  REAL(KIND=8), INTENT(IN) :: s
334  REAL(KIND=8) :: y
335  INTEGER :: i,j
336  corr_ij=0.d0
337  DO i=1,n_unres
338  DO j=1,n_unres
339  CALL splint(xa,ya(i,j,:),y2(i,j,:),nspl,s,y)
340  corr_ij(i,j)=y
341  END DO
342  END DO
real(kind=8) function corrmod::fs ( real(kind=8), intent(in)  s,
real(kind=8), dimension(5), intent(in)  p 
)
private

Exponential fit function.

Parameters
stime \(s\) at which the function is evaluated
pvector holding the coefficients of the fit function

Definition at line 388 of file corrmod.f90.

388  REAL(KIND=8), INTENT(IN) :: s
389  REAL(KIND=8), DIMENSION(5), INTENT(IN) :: p
390  REAL(KIND=8) :: fs
391  fs=p(1)*exp(-s/p(2))*cos(p(3)*s+p(4))
392  RETURN
subroutine, public corrmod::init_corr ( )

Subroutine to initialise the computation of the correlation.

Definition at line 46 of file corrmod.f90.

46  INTEGER :: allocstat,i,j,k,nf
47  REAL(KIND=8), DIMENSION(5) :: dumb
48  LOGICAL :: ex
49 
50  ! Selection of the loading mode
51  SELECT CASE (load_mode)
52  CASE ('defi')
53  corrcomp => corrcomp_from_def
54  CASE ('spli')
55  INQUIRE(file='corrspline.def',exist=ex)
56  IF (.not.ex) stop "*** File corrspline.def not found ! ***"
57  OPEN(20,file='corrspline.def',status='old')
58  READ(20,*) nf,nspl
59  IF (nf /= n_unres) stop "*** Dimension in files corrspline.def and sf.nml do not correspond ! ***"
60  ALLOCATE(xa(nspl), ya(n_unres,n_unres,nspl), y2(n_unres,n_unres,nspl), stat=allocstat)
61  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
62  READ(20,*) xa
63  maxint=xa(nspl)/2
64  DO k=1,n_unres*n_unres
65  READ(20,*) i,j
66  READ(20,*) ya(i,j,:)
67  READ(20,*) y2(i,j,:)
68  ENDDO
69  CLOSE(20)
70  corrcomp => corrcomp_from_spline
71  klo=1
72  khi=nspl
73  CASE ('expo')
74  INQUIRE(file='correxpo.def',exist=ex)
75  IF (.not.ex) stop "*** File correxpo.def not found ! ***"
76  OPEN(20,file='correxpo.def',status='old')
77  READ(20,*) nf,maxint
78  IF (nf /= n_unres) stop "*** Dimension in files correxpo.def and sf.nml do not correspond ! ***"
79  ALLOCATE(ya(n_unres,n_unres,5), stat=allocstat)
80  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
81  DO k=1,n_unres*n_unres
82  READ(20,*) i,j,dumb
83  ya(i,j,:)=dumb
84  ENDDO
85  CLOSE(20)
86  corrcomp => corrcomp_from_fit
87  CASE DEFAULT
88  stop '*** LOAD_MODE variable not properly defined in corrmod.nml ***'
89  END SELECT
90 
91  ALLOCATE(mean(n_unres),mean_full(0:ndim), stat=allocstat)
92  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
93 
94  ALLOCATE(inv_corr_i(n_unres,n_unres), stat=allocstat)
95  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
96 
97  ALLOCATE(corr_i(n_unres,n_unres), stat=allocstat)
98  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
99 
100  ALLOCATE(corr_ij(n_unres,n_unres), stat=allocstat)
101  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
102 
103  ALLOCATE(corr_i_full(ndim,ndim), stat=allocstat)
104  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
105 
106  ALLOCATE(inv_corr_i_full(ndim,ndim), stat=allocstat)
107  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
108 
109  corr_ij=0.d0
110 
111  CALL corrcomp(0.d0)
112  corr_i=corr_ij
113  inv_corr_i=invmat(corr_i)
114 
115  corr_i_full=0.d0
116  DO i=1,n_unres
117  DO j=1,n_unres
118  corr_i_full(ind(i),ind(j))=corr_i(i,j)
119  ENDDO
120  ENDDO
121 
122  inv_corr_i_full=0.d0
123  DO i=1,n_unres
124  DO j=1,n_unres
125  inv_corr_i_full(ind(i),ind(j))=inv_corr_i(i,j)
126  ENDDO
127  ENDDO
128 
129  mean=0.d0
130  INQUIRE(file='mean.def',exist=ex)
131  IF (ex) THEN
132  OPEN(20,file='mean.def',status='old')
133  READ(20,*) mean
134  CLOSE(20)
135  ENDIF
136 
137  mean_full=0.d0
138  DO i=1,n_unres
139  mean_full(ind(i))=mean(i)
140  ENDDO
141 
Statistics accumulators.
Definition: stat.f90:14
subroutine corrmod::splint ( real(kind=8), dimension(n), intent(in)  xa,
real(kind=8), dimension(n), intent(in)  ya,
real(kind=8), dimension(n), intent(in)  y2a,
integer, intent(in)  n,
real(kind=8), intent(in)  x,
real(kind=8), intent(out)  y 
)
private

Routine to compute the spline representation parameters.

Definition at line 347 of file corrmod.f90.

347  INTEGER, INTENT(IN) :: n
348  REAL(KIND=8), INTENT(IN), DIMENSION(n) :: xa,y2a,ya
349  REAL(KIND=8), INTENT(IN) :: x
350  REAL(KIND=8), INTENT(OUT) :: y
351  INTEGER :: k
352  REAL(KIND=8) :: a,b,h
353  if ((khi-klo.gt.1).or.(xa(klo).gt.x).or.(xa(khi).lt.x)) then
354  if ((khi-klo.eq.1).and.(xa(klo).lt.x)) then
355  khi=klo
356  DO WHILE (xa(khi).lt.x)
357  khi=khi+1
358  END DO
359  klo=khi-1
360  else
361  khi=n
362  klo=1
363  DO WHILE (khi-klo.gt.1)
364  k=(khi+klo)/2
365  if(xa(k).gt.x)then
366  khi=k
367  else
368  klo=k
369  endif
370  END DO
371  end if
372  ! print*, "search",x,khi-klo,xa(klo),xa(khi)
373  ! else
374  ! print*, "ok",x,khi-klo,xa(klo),xa(khi)
375  endif
376  h=xa(khi)-xa(klo)
377  if (h.eq.0.) stop 'bad xa input in splint'
378  a=(xa(khi)-x)/h
379  b=(x-xa(klo))/h
380  y=a*ya(klo)+b*ya(khi)+((a**3-a)*y2a(klo)+(b**3-b)*y2a(khi))*(h**2)/6.
381  return

Variable Documentation

real(kind=8), dimension(:,:), allocatable, public corrmod::corr_i

Covariance matrix of the unresolved variables (reduced version)

Definition at line 30 of file corrmod.f90.

30  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE, PUBLIC :: corr_i !< Covariance matrix of the unresolved variables (reduced version)
real(kind=8), dimension(:,:), allocatable, public corrmod::corr_i_full

Covariance matrix of the unresolved variables (full version)

Definition at line 28 of file corrmod.f90.

28  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE, PUBLIC :: corr_i_full !< Covariance matrix of the unresolved variables (full version)
real(kind=8), dimension(:,:), allocatable, public corrmod::corr_ij

Matrix holding the correlation matrix at a given time.

Definition at line 32 of file corrmod.f90.

32  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE, PUBLIC :: corr_ij !< Matrix holding the correlation matrix at a given time
procedure(corrcomp_from_spline), pointer, public corrmod::corrcomp

Pointer to the correlation computation routine.

Definition at line 41 of file corrmod.f90.

41  PROCEDURE(corrcomp_from_spline), POINTER, PUBLIC :: corrcomp
real(kind=8), dimension(:,:), allocatable, public corrmod::inv_corr_i

Inverse of the covariance matrix of the unresolved variables (reduced version)

Definition at line 31 of file corrmod.f90.

31  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE, PUBLIC :: inv_corr_i !< Inverse of the covariance matrix of the unresolved variables (reduced version)
real(kind=8), dimension(:,:), allocatable, public corrmod::inv_corr_i_full

Inverse of the covariance matrix of the unresolved variables (full version)

Definition at line 29 of file corrmod.f90.

29  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE, PUBLIC :: inv_corr_i_full !< Inverse of the covariance matrix of the unresolved variables (full version)
integer corrmod::khi
private

Definition at line 38 of file corrmod.f90.

integer corrmod::klo
private

Definition at line 38 of file corrmod.f90.

real(kind=8), dimension(:), allocatable, public corrmod::mean

Vector holding the mean of the unresolved dynamics (reduced version)

Definition at line 26 of file corrmod.f90.

26  REAL(KIND=8), DIMENSION(:), ALLOCATABLE, PUBLIC :: mean !< Vector holding the mean of the unresolved dynamics (reduced version)
real(kind=8), dimension(:), allocatable, public corrmod::mean_full

Vector holding the mean of the unresolved dynamics (full version)

Definition at line 27 of file corrmod.f90.

27  REAL(KIND=8), DIMENSION(:), ALLOCATABLE, PUBLIC :: mean_full !< Vector holding the mean of the unresolved dynamics (full version)
integer corrmod::nspl
private

Integers needed by the spline representation of the correlation.

Definition at line 38 of file corrmod.f90.

38  INTEGER :: nspl,klo,khi
real(kind=8), dimension(:), allocatable corrmod::xa
private

Vector holding coefficient of the spline and exponential correlation representation.

Definition at line 35 of file corrmod.f90.

35  REAL(KIND=8), DIMENSION(:), ALLOCATABLE :: xa !< Vector holding coefficient of the spline and exponential correlation representation
real(kind=8), dimension(:,:,:), allocatable corrmod::y2
private

Vector holding coefficient of the spline and exponential correlation representation.

Definition at line 33 of file corrmod.f90.

33  REAL(KIND=8), DIMENSION(:,:,:), ALLOCATABLE :: y2 !< Vector holding coefficient of the spline and exponential correlation representation
real(kind=8), dimension(:,:,:), allocatable corrmod::ya
private

Vector holding coefficient of the spline and exponential correlation representation.

Definition at line 34 of file corrmod.f90.

34  REAL(KIND=8), DIMENSION(:,:,:), ALLOCATABLE :: ya !< Vector holding coefficient of the spline and exponential correlation representation