EFTCAMB  Reference documentation for version 3.0
07p1_Pure_EFT_std.f90
Go to the documentation of this file.
1 !----------------------------------------------------------------------------------------
2 !
3 ! This file is part of EFTCAMB.
4 !
5 ! Copyright (C) 2013-2016 by the EFTCAMB authors
6 !
7 ! The EFTCAMB code is free software;
8 ! You can use it, redistribute it, and/or modify it under the terms
9 ! of the GNU General Public License as published by the Free Software Foundation;
10 ! either version 3 of the License, or (at your option) any later version.
11 ! The full text of the license can be found in the file eftcamb/LICENSE at
12 ! the top level of the EFTCAMB distribution.
13 !
14 !----------------------------------------------------------------------------------------
15 
19 
20 
21 !----------------------------------------------------------------------------------------
24 
26 
28 
29  use precision
30  use inifile
31  use amlutils
32  use eftcamb_cache
33  use eft_def
45 
46  implicit none
47 
48  private
49 
50  public eftcamb_std_pure_eft
51 
52  !----------------------------------------------------------------------------------------
54  type, extends ( eftcamb_designer_model ) :: eftcamb_std_pure_eft
55 
56  ! the pure EFT functions model selection flags:
57  integer :: pureeftmodelomega
58  integer :: eftwde
59  integer :: pureeftmodelgamma1
60  integer :: pureeftmodelgamma2
61  integer :: pureeftmodelgamma3
62  integer :: pureeftmodelgamma4
63  integer :: pureeftmodelgamma5
64  integer :: pureeftmodelgamma6
65 
66  ! selection flag for Horndeski:
67  logical :: pureefthorndeski
68 
69  ! the pure EFT functions:
70  class( parametrized_function_1d ), allocatable :: pureeftomega
71  class( parametrized_function_1d ), allocatable :: pureeftwde
72  class( parametrized_function_1d ), allocatable :: pureeftgamma1
73  class( parametrized_function_1d ), allocatable :: pureeftgamma2
74  class( parametrized_function_1d ), allocatable :: pureeftgamma3
75  class( parametrized_function_1d ), allocatable :: pureeftgamma4
76  class( parametrized_function_1d ), allocatable :: pureeftgamma5
77  class( parametrized_function_1d ), allocatable :: pureeftgamma6
78 
79  contains
80 
81  ! initialization of the model:
82  procedure :: read_model_selection => eftcambpureeftstdreadmodelselectionfromfile
83  procedure :: allocate_model_selection => eftcambpureeftstdallocatemodelselection
84  procedure :: init_model_parameters => eftcambpureeftstdinitmodelparameters
85  procedure :: init_model_parameters_from_file => eftcambpureeftstdinitmodelparametersfromfile
86  ! utility functions:
87  procedure :: compute_param_number => eftcambpureeftstdcomputeparametersnumber
88  procedure :: feedback => eftcambpureeftstdfeedback
89  procedure :: parameter_names => eftcambpureeftstdparameternames
90  procedure :: parameter_names_latex => eftcambpureeftstdparameternameslatex
91  procedure :: parameter_values => eftcambpureeftstdparametervalues
92  ! CAMB related procedures:
93  procedure :: compute_background_eft_functions => eftcambpureeftstdbackgroundeftfunctions
94  procedure :: compute_secondorder_eft_functions => eftcambpureeftstdsecondordereftfunctions
95  procedure :: compute_dtauda => eftcambpureeftstdcomputedtauda
96  procedure :: compute_adotoa => eftcambpureeftstdcomputeadotoa
97  procedure :: compute_h_derivs => eftcambpureeftstdcomputehubbleder
98  ! stability procedures:
99  procedure :: additional_model_stability => eftcambpureeftstdadditionalmodelstability
100 
101  end type eftcamb_std_pure_eft
102 
103  ! ---------------------------------------------------------------------------------------------
104 
105 contains
106 
107  ! ---------------------------------------------------------------------------------------------
109  subroutine eftcambpureeftstdreadmodelselectionfromfile( self, Ini )
110 
111  implicit none
112 
113  class(eftcamb_std_pure_eft) :: self
114  type(tinifile) :: ini
115 
116  ! read model selection flags:
117  self%PureEFTmodelOmega = ini_read_int_file( ini, 'PureEFTmodelOmega' , 0 )
118  self%EFTwDE = ini_read_int_file( ini, 'EFTwDE' , 0 )
119  self%PureEFTmodelGamma1 = ini_read_int_file( ini, 'PureEFTmodelGamma1' , 0 )
120  self%PureEFTmodelGamma2 = ini_read_int_file( ini, 'PureEFTmodelGamma2' , 0 )
121  self%PureEFTmodelGamma3 = ini_read_int_file( ini, 'PureEFTmodelGamma3' , 0 )
122  self%PureEFTmodelGamma4 = ini_read_int_file( ini, 'PureEFTmodelGamma4' , 0 )
123  self%PureEFTmodelGamma5 = ini_read_int_file( ini, 'PureEFTmodelGamma5' , 0 )
124  self%PureEFTmodelGamma6 = ini_read_int_file( ini, 'PureEFTmodelGamma6' , 0 )
125  ! read the Horndeski flag:
126  self%PureEFTHorndeski = ini_read_logical_file( ini, 'PureEFTHorndeski' , .false. )
127 
128  end subroutine eftcambpureeftstdreadmodelselectionfromfile
129 
130  ! ---------------------------------------------------------------------------------------------
132  subroutine eftcambpureeftstdallocatemodelselection( self )
133 
134  implicit none
135 
136  class(eftcamb_std_pure_eft) :: self
137 
138  ! allocate Omega:
139  if ( allocated(self%PureEFTOmega) ) deallocate(self%PureEFTOmega)
140  select case ( self%PureEFTmodelOmega )
141  case(0)
142  allocate( zero_parametrization_1d::self%PureEFTOmega )
143  case(1)
144  allocate( constant_parametrization_1d::self%PureEFTOmega )
145  case(2)
146  allocate( linear_parametrization_1d::self%PureEFTOmega )
147  case(3)
148  allocate( power_law_parametrization_1d::self%PureEFTOmega )
149  call self%PureEFTOmega%set_param_names( ['EFTOmega0 ', 'EFTOmegaExp'], ['\Omega_0^{\rm EFT}', 'n^{\rm EFT} '] )
150  case(4)
151  allocate( exponential_parametrization_1d::self%PureEFTOmega )
152  call self%PureEFTOmega%set_param_names( ['EFTOmega0 ', 'EFTOmegaExp'], ['\Omega_0^{\rm EFT}', 'n^{\rm EFT} '] )
153  case default
154  write(*,'(a,I3)') 'No model corresponding to PureEFTmodelOmega =', self%PureEFTmodelOmega
155  write(*,'(a)') 'Please select an appropriate model.'
156  end select
157  ! allocate wDE:
158  if ( allocated(self%PureEFTwDE) ) deallocate(self%PureEFTwDE)
159  select case ( self%EFTwDE )
160  case(0)
161  allocate( wde_lcdm_parametrization_1d::self%PureEFTwDE )
162  case(1)
163  allocate( constant_parametrization_1d::self%PureEFTwDE )
164  case(2)
165  allocate( cpl_parametrization_1d::self%PureEFTwDE )
166  call self%PureEFTwDE%set_param_names( ['EFTw0', 'EFTwa'], ['w_0', 'w_a'] )
167  case(3)
168  allocate( jbp_parametrization_1d::self%PureEFTwDE )
169  call self%PureEFTwDE%set_param_names( ['EFTw0', 'EFTwa', 'EFTwn'], [ 'w_0', 'w_a', 'n ' ] )
170  case(4)
171  allocate( turning_point_parametrization_1d::self%PureEFTwDE )
172  call self%PureEFTwDE%set_param_names( ['EFTw0 ', 'EFTwa ', 'EFTwat'], ['w_0', 'w_a', 'a_t'] )
173  case(5)
174  allocate( taylor_parametrization_1d::self%PureEFTwDE )
175  call self%PureEFTwDE%set_param_names( ['EFTw0', 'EFTwa', 'EFTw2', 'EFTw3'], ['w_0', 'w_a', 'w_2', 'w_3'] )
176  case default
177  write(*,'(a,I3)') 'No model corresponding to EFTwDE =', self%EFTwDE
178  write(*,'(a)') 'Please select an appropriate model.'
179  end select
180  ! allocate Gamma1:
181  if ( allocated(self%PureEFTGamma1) ) deallocate(self%PureEFTGamma1)
182  select case ( self%PureEFTmodelGamma1 )
183  case(0)
184  allocate( zero_parametrization_1d::self%PureEFTGamma1 )
185  case(1)
186  allocate( constant_parametrization_1d::self%PureEFTGamma1 )
187  case(2)
188  allocate( linear_parametrization_1d::self%PureEFTGamma1 )
189  case(3)
190  allocate( power_law_parametrization_1d::self%PureEFTGamma1 )
191  call self%PureEFTGamma1%set_param_names( ['EFTGamma10 ', 'EFTGamma1Exp'], ['\gamma_0^{(1) {\rm EFT}} ', '\gamma_{\rm exp}^{(1) {\rm EFT}}'] )
192  case(4)
193  allocate( exponential_parametrization_1d::self%PureEFTGamma1 )
194  call self%PureEFTGamma1%set_param_names( ['EFTGamma10 ', 'EFTGamma1Exp'], ['\gamma_0^{(1) {\rm EFT}} ', '\gamma_{\rm exp}^{(1) {\rm EFT}}'] )
195  case default
196  write(*,'(a,I3)') 'No model corresponding to PureEFTmodelGamma1 =', self%PureEFTmodelGamma1
197  write(*,'(a)') 'Please select an appropriate model.'
198  end select
199  ! allocate Gamma2:
200  if ( allocated(self%PureEFTGamma2) ) deallocate(self%PureEFTGamma2)
201  select case ( self%PureEFTmodelGamma2 )
202  case(0)
203  allocate( zero_parametrization_1d::self%PureEFTGamma2 )
204  case(1)
205  allocate( constant_parametrization_1d::self%PureEFTGamma2 )
206  case(2)
207  allocate( linear_parametrization_1d::self%PureEFTGamma2 )
208  case(3)
209  allocate( power_law_parametrization_1d::self%PureEFTGamma2 )
210  call self%PureEFTGamma2%set_param_names( ['EFTGamma20 ', 'EFTGamma2Exp'], ['\gamma_0^{(2) {\rm EFT}} ', '\gamma_{\rm exp}^{(2) {\rm EFT}}'] )
211  case(4)
212  allocate( exponential_parametrization_1d::self%PureEFTGamma2 )
213  call self%PureEFTGamma2%set_param_names( ['EFTGamma20 ', 'EFTGamma2Exp'], ['\gamma_0^{(2) {\rm EFT}} ', '\gamma_{\rm exp}^{(2) {\rm EFT}}'] )
214  case default
215  write(*,'(a,I3)') 'No model corresponding to PureEFTmodelGamma2 =', self%PureEFTmodelGamma2
216  write(*,'(a)') 'Please select an appropriate model.'
217  end select
218  ! allocate Gamma3:
219  if ( allocated(self%PureEFTGamma3) ) deallocate(self%PureEFTGamma3)
220  select case ( self%PureEFTmodelGamma3 )
221  case(0)
222  allocate( zero_parametrization_1d::self%PureEFTGamma3 )
223  case(1)
224  allocate( constant_parametrization_1d::self%PureEFTGamma3 )
225  case(2)
226  allocate( linear_parametrization_1d::self%PureEFTGamma3 )
227  case(3)
228  allocate( power_law_parametrization_1d::self%PureEFTGamma3 )
229  call self%PureEFTGamma3%set_param_names( ['EFTGamma30 ', 'EFTGamma3Exp'], ['\gamma_0^{(3) {\rm EFT}} ', '\gamma_{\rm exp}^{(3) {\rm EFT}}'] )
230  case(4)
231  allocate( exponential_parametrization_1d::self%PureEFTGamma3 )
232  call self%PureEFTGamma3%set_param_names( ['EFTGamma30 ', 'EFTGamma3Exp'], ['\gamma_0^{(3) {\rm EFT}} ', '\gamma_{\rm exp}^{(3) {\rm EFT}}'] )
233  case default
234  write(*,'(a,I3)') 'No model corresponding to PureEFTmodelGamma3 =', self%PureEFTmodelGamma3
235  write(*,'(a)') 'Please select an appropriate model.'
236  end select
237  ! allocate the other functions only if not Horndeski:
238  if ( .not. self%PureEFTHorndeski ) then
239  ! allocate Gamma4:
240  if ( allocated(self%PureEFTGamma4) ) deallocate(self%PureEFTGamma4)
241  select case ( self%PureEFTmodelGamma4 )
242  case(0)
243  allocate( zero_parametrization_1d::self%PureEFTGamma4 )
244  case(1)
245  allocate( constant_parametrization_1d::self%PureEFTGamma4 )
246  case(2)
247  allocate( linear_parametrization_1d::self%PureEFTGamma4 )
248  case(3)
249  allocate( power_law_parametrization_1d::self%PureEFTGamma4 )
250  call self%PureEFTGamma4%set_param_names( ['EFTGamma40 ', 'EFTGamma4Exp'], ['\gamma_0^{(4) {\rm EFT}} ', '\gamma_{\rm exp}^{(4) {\rm EFT}}'] )
251  case(4)
252  allocate( exponential_parametrization_1d::self%PureEFTGamma4 )
253  call self%PureEFTGamma4%set_param_names( ['EFTGamma40 ', 'EFTGamma4Exp'], ['\gamma_0^{(4) {\rm EFT}} ', '\gamma_{\rm exp}^{(4) {\rm EFT}}'] )
254  case default
255  write(*,'(a,I3)') 'No model corresponding to PureEFTmodelGamma4 =', self%PureEFTmodelGamma4
256  write(*,'(a)') 'Please select an appropriate model.'
257  end select
258  ! allocate Gamma5:
259  if ( allocated(self%PureEFTGamma5) ) deallocate(self%PureEFTGamma5)
260  select case ( self%PureEFTmodelGamma5 )
261  case(0)
262  allocate( zero_parametrization_1d::self%PureEFTGamma5 )
263  case(1)
264  allocate( constant_parametrization_1d::self%PureEFTGamma5 )
265  case(2)
266  allocate( linear_parametrization_1d::self%PureEFTGamma5 )
267  case(3)
268  allocate( power_law_parametrization_1d::self%PureEFTGamma5 )
269  call self%PureEFTGamma5%set_param_names( ['EFTGamma50 ', 'EFTGamma5Exp'], ['\gamma_0^{(5) {\rm EFT}} ', '\gamma_{\rm exp}^{(5) {\rm EFT}}'] )
270  case(4)
271  allocate( exponential_parametrization_1d::self%PureEFTGamma5 )
272  call self%PureEFTGamma5%set_param_names( ['EFTGamma50 ', 'EFTGamma5Exp'], ['\gamma_0^{(5) {\rm EFT}} ', '\gamma_{\rm exp}^{(5) {\rm EFT}}'] )
273  case default
274  write(*,'(a,I3)') 'No model corresponding to PureEFTmodelGamma5 =', self%PureEFTmodelGamma5
275  write(*,'(a)') 'Please select an appropriate model.'
276  end select
277  ! allocate Gamma6:
278  if ( allocated(self%PureEFTGamma6) ) deallocate(self%PureEFTGamma6)
279  select case ( self%PureEFTmodelGamma6 )
280  case(0)
281  allocate( zero_parametrization_1d::self%PureEFTGamma6 )
282  case(1)
283  allocate( constant_parametrization_1d::self%PureEFTGamma6 )
284  case(2)
285  allocate( linear_parametrization_1d::self%PureEFTGamma6 )
286  case(3)
287  allocate( power_law_parametrization_1d::self%PureEFTGamma6 )
288  call self%PureEFTGamma6%set_param_names( ['EFTGamma60 ', 'EFTGamma6Exp'], ['\gamma_0^{(6) {\rm EFT}} ', '\gamma_{\rm exp}^{(6) {\rm EFT}}'] )
289  case(4)
290  allocate( exponential_parametrization_1d::self%PureEFTGamma6 )
291  call self%PureEFTGamma6%set_param_names( ['EFTGamma60 ', 'EFTGamma6Exp'], ['\gamma_0^{(6) {\rm EFT}} ', '\gamma_{\rm exp}^{(6) {\rm EFT}}'] )
292  case default
293  write(*,'(a,I3)') 'No model corresponding to PureEFTmodelGamma6 =', self%PureEFTmodelGamma6
294  write(*,'(a)') 'Please select an appropriate model.'
295  end select
296  end if
297 
298  ! initialize the names:
299  call self%PureEFTOmega%set_name ( 'EFTOmega' , '\Omega' )
300  call self%PureEFTwDE%set_name ( 'EFTw' , 'w' )
301  call self%PureEFTGamma1%set_name( 'EFTGamma1', '\gamma^{(1)}' )
302  call self%PureEFTGamma2%set_name( 'EFTGamma2', '\gamma^{(2)}' )
303  call self%PureEFTGamma3%set_name( 'EFTGamma3', '\gamma^{(3)}' )
304 
305  if ( .not. self%PureEFTHorndeski ) then
306  call self%PureEFTGamma4%set_name( 'EFTGamma4', '\gamma^{(4)}' )
307  call self%PureEFTGamma5%set_name( 'EFTGamma5', '\gamma^{(5)}' )
308  call self%PureEFTGamma6%set_name( 'EFTGamma6', '\gamma^{(6)}' )
309  end if
310 
311  end subroutine eftcambpureeftstdallocatemodelselection
312 
313  ! ---------------------------------------------------------------------------------------------
315  subroutine eftcambpureeftstdinitmodelparameters( self, array )
316 
317  implicit none
318 
319  class(eftcamb_std_pure_eft) :: self
320  real(dl), dimension(self%parameter_number), intent(in) :: array
321 
322  real(dl), allocatable, dimension(:) :: temp
323  integer :: num_params_function, num_params_temp, i
324 
325  num_params_temp = 1
326 
327  ! first elements are Omega parameters:
328  num_params_function = self%PureEFTOmega%parameter_number
329  allocate( temp(num_params_function) )
330  do i = 1, num_params_function
331  temp(i) = array(num_params_temp)
332  num_params_temp = num_params_temp +1
333  end do
334  call self%PureEFTOmega%init_parameters(temp)
335  deallocate( temp )
336  ! then w_DE parameters:
337  num_params_function = self%PureEFTwDE%parameter_number
338  allocate( temp(num_params_function) )
339  do i = 1, num_params_function
340  temp(i) = array(num_params_temp)
341  num_params_temp = num_params_temp +1
342  end do
343  call self%PureEFTwDE%init_parameters(temp)
344  deallocate(temp)
345  ! then gamma1:
346  num_params_function = self%PureEFTGamma1%parameter_number
347  allocate( temp(num_params_function) )
348  do i = 1, num_params_function
349  temp(i) = array(num_params_temp)
350  num_params_temp = num_params_temp +1
351  end do
352  call self%PureEFTGamma1%init_parameters(temp)
353  deallocate(temp)
354  ! then gamma2:
355  num_params_function = self%PureEFTGamma2%parameter_number
356  allocate( temp(num_params_function) )
357  do i = 1, num_params_function
358  temp(i) = array(num_params_temp)
359  num_params_temp = num_params_temp +1
360  end do
361  call self%PureEFTGamma2%init_parameters(temp)
362  deallocate(temp)
363  ! then gamma3:
364  num_params_function = self%PureEFTGamma3%parameter_number
365  allocate( temp(num_params_function) )
366  do i = 1, num_params_function
367  temp(i) = array(num_params_temp)
368  num_params_temp = num_params_temp +1
369  end do
370  call self%PureEFTGamma3%init_parameters(temp)
371  deallocate(temp)
372 
373  ! then beyond Horndeski parameters:
374  if ( .not. self%PureEFTHorndeski ) then
375 
376  ! gamma4:
377  num_params_function = self%PureEFTGamma4%parameter_number
378  allocate( temp(num_params_function) )
379  do i = 1, num_params_function
380  temp(i) = array(num_params_temp)
381  num_params_temp = num_params_temp +1
382  end do
383  call self%PureEFTGamma4%init_parameters(temp)
384  deallocate(temp)
385  ! gamma5:
386  num_params_function = self%PureEFTGamma5%parameter_number
387  allocate( temp(num_params_function) )
388  do i = 1, num_params_function
389  temp(i) = array(num_params_temp)
390  num_params_temp = num_params_temp +1
391  end do
392  call self%PureEFTGamma5%init_parameters(temp)
393  deallocate(temp)
394  ! gamma6:
395  num_params_function = self%PureEFTGamma6%parameter_number
396  allocate( temp(num_params_function) )
397  do i = 1, num_params_function
398  temp(i) = array(num_params_temp)
399  num_params_temp = num_params_temp +1
400  end do
401  call self%PureEFTGamma6%init_parameters(temp)
402  deallocate(temp)
403 
404  end if
405 
406  ! now check the length of the parameters:
407  if ( num_params_temp-1 /= self%parameter_number ) then
408  write(*,*) 'In EFTCAMBPureEFTstdInitModelParameters:'
409  write(*,*) 'Length of num_params_temp and self%parameter_number do not coincide.'
410  write(*,*) 'num_params_temp:', num_params_temp-1
411  write(*,*) 'self%parameter_number:', self%parameter_number
412  call mpistop('EFTCAMB error')
413  end if
414 
415  end subroutine eftcambpureeftstdinitmodelparameters
416 
417  ! ---------------------------------------------------------------------------------------------
419  subroutine eftcambpureeftstdinitmodelparametersfromfile( self, Ini )
420 
421  implicit none
422 
423  class(eftcamb_std_pure_eft) :: self
424  type(tinifile) :: ini
425 
426  call self%PureEFTOmega%init_from_file ( ini )
427  call self%PureEFTwDE%init_from_file ( ini )
428  call self%PureEFTGamma1%init_from_file( ini )
429  call self%PureEFTGamma2%init_from_file( ini )
430  call self%PureEFTGamma3%init_from_file( ini )
431 
432  if ( .not. self%PureEFTHorndeski ) then
433  call self%PureEFTGamma4%init_from_file( ini )
434  call self%PureEFTGamma5%init_from_file( ini )
435  call self%PureEFTGamma6%init_from_file( ini )
436  end if
437 
438  end subroutine eftcambpureeftstdinitmodelparametersfromfile
439 
440  ! ---------------------------------------------------------------------------------------------
442  subroutine eftcambpureeftstdcomputeparametersnumber( self )
443 
444  implicit none
445 
446  class(eftcamb_std_pure_eft) :: self
447 
448  self%parameter_number = 0
449  self%parameter_number = self%parameter_number +self%PureEFTOmega%parameter_number
450  self%parameter_number = self%parameter_number +self%PureEFTwDE%parameter_number
451  self%parameter_number = self%parameter_number +self%PureEFTGamma1%parameter_number
452  self%parameter_number = self%parameter_number +self%PureEFTGamma2%parameter_number
453  self%parameter_number = self%parameter_number +self%PureEFTGamma3%parameter_number
454  if ( .not. self%PureEFTHorndeski ) then
455  self%parameter_number = self%parameter_number +self%PureEFTGamma4%parameter_number
456  self%parameter_number = self%parameter_number +self%PureEFTGamma5%parameter_number
457  self%parameter_number = self%parameter_number +self%PureEFTGamma6%parameter_number
458  end if
459 
460  end subroutine eftcambpureeftstdcomputeparametersnumber
461 
462  ! ---------------------------------------------------------------------------------------------
464  subroutine eftcambpureeftstdfeedback( self, print_params )
465 
466  implicit none
467 
468  class(eftcamb_std_pure_eft) :: self
469  logical, optional :: print_params
471 
472  ! print general model informations:
473  write(*,*)
474  write(*,'(a,a)') ' Model = ', self%name
475  if ( self%PureEFTHorndeski ) then
476  write(*,"(a)") ' Pure EFT Horndeski'
477  end if
478  write(*,'(a,I3)') ' Number of params =' , self%parameter_number
479 
480  ! print model functions informations:
481  write(*,*)
482  if ( self%PureEFTmodelOmega /= 0 ) write(*,'(a,I3)') ' PureEFTmodelOmega =', self%PureEFTmodelOmega
483  if ( self%EFTwDE /= 0 ) write(*,'(a,I3)') ' EFTwDE =', self%EFTwDE
484  if ( self%PureEFTmodelGamma1 /= 0 ) write(*,'(a,I3)') ' PureEFTmodelGamma1 =', self%PureEFTmodelGamma1
485  if ( self%PureEFTmodelGamma2 /= 0 ) write(*,'(a,I3)') ' PureEFTmodelGamma2 =', self%PureEFTmodelGamma2
486  if ( self%PureEFTmodelGamma3 /= 0 ) write(*,'(a,I3)') ' PureEFTmodelGamma3 =', self%PureEFTmodelGamma3
487  if ( .not. self%PureEFTHorndeski ) then
488  if ( self%PureEFTmodelGamma4 /= 0 ) write(*,'(a,I3)') ' PureEFTmodelGamma4 =', self%PureEFTmodelGamma4
489  if ( self%PureEFTmodelGamma5 /= 0 ) write(*,'(a,I3)') ' PureEFTmodelGamma5 =', self%PureEFTmodelGamma5
490  if ( self%PureEFTmodelGamma6 /= 0 ) write(*,'(a,I3)') ' PureEFTmodelGamma6 =', self%PureEFTmodelGamma6
491  end if
492 
493  write(*,*)
494  ! print functions informations:
495  call self%PureEFTOmega%feedback ( print_params )
496  call self%PureEFTwDE%feedback ( print_params )
497  call self%PureEFTGamma1%feedback ( print_params )
498  call self%PureEFTGamma2%feedback ( print_params )
499  call self%PureEFTGamma3%feedback ( print_params )
500  if ( .not. self%PureEFTHorndeski ) then
501  call self%PureEFTGamma4%feedback( print_params )
502  call self%PureEFTGamma5%feedback( print_params )
503  call self%PureEFTGamma6%feedback( print_params )
504  end if
505 
506  end subroutine eftcambpureeftstdfeedback
507 
508  ! ---------------------------------------------------------------------------------------------
510  subroutine eftcambpureeftstdparameternames( self, i, name )
511 
512  implicit none
513 
514  class(eftcamb_std_pure_eft) :: self
515  integer , intent(in) :: i
516  character(*), intent(out) :: name
517 
518  integer :: nomega, nw, n1, n2, n3, n4, n5, n6
519  integer :: j
520 
521  ! compute the incremental number of parameters:
522  nomega = self%PureEFTOmega%parameter_number
523  nw = nomega + self%PureEFTwDE%parameter_number
524  n1 = nw + self%PureEFTGamma1%parameter_number
525  n2 = n1 + self%PureEFTGamma2%parameter_number
526  n3 = n2 + self%PureEFTGamma3%parameter_number
527  if ( .not. self%PureEFTHorndeski ) then
528  n4 = n3 + self%PureEFTGamma4%parameter_number
529  n5 = n4 + self%PureEFTGamma5%parameter_number
530  n6 = n5 + self%PureEFTGamma6%parameter_number
531  end if
532 
533  ! check validity of input:
534  if ( i > self%parameter_number .or. i <= 0 ) then
535  write(*,'(a,I3)') 'No parameter corresponding to: ', i
536  write(*,'(a,I3)') 'Total number of parameters is: ', self%parameter_number
537  call mpistop('EFTCAMB error')
538 
539  ! parameter from Omega function
540  else if ( i <= nomega ) then
541  do j = 1, self%PureEFTOmega%parameter_number
542  if ( i == j ) call self%PureEFTOmega%parameter_names( j, name )
543  end do
544  return
545 
546  ! parameter from wDE function
547  else if ( i <= nw ) then
548  do j = 1, self%PureEFTwDE%parameter_number
549  if ( i-nomega == j ) call self%PureEFTwDE%parameter_names( j, name )
550  end do
551  return
552 
553  ! parameter from Gamma1 function
554  else if ( i <= n1) then
555  do j = 1, self%PureEFTGamma1%parameter_number
556  if ( i-nw == j ) call self%PureEFTGamma1%parameter_names( j, name )
557  end do
558  return
559 
560  !parameter from Gamma2 function
561  else if ( i <= n2 ) then
562  do j = 1, self%PureEFTGamma2%parameter_number
563  if ( i-n1 == j ) call self%PureEFTGamma2%parameter_names( j, name )
564  end do
565  return
566 
567  !parameter from Gamma3 function
568  else if ( i <= n3 ) then
569 
570  do j = 1, self%PureEFTGamma3%parameter_number
571  if ( i-n2 == j ) call self%PureEFTGamma3%parameter_names( j, name )
572  end do
573  return
574 
575  !parameter from Gamma4 function
576  else if ( .not. self%PureEFTHorndeski .and. i <= n4 ) then
577 
578  do j = 1, self%PureEFTGamma4%parameter_number
579  if ( i-n3 == j ) call self%PureEFTGamma4%parameter_names( j, name )
580  end do
581  return
582 
583  !parameter from Gamma5 function
584  else if ( .not. self%PureEFTHorndeski .and. i <= n5 ) then
585 
586  do j = 1, self%PureEFTGamma5%parameter_number
587  if ( i-n4 == j ) call self%PureEFTGamma5%parameter_names( j, name )
588  end do
589  return
590 
591  !parameter from Gamma6 function
592  else if ( .not. self%PureEFTHorndeski .and. i <= n6 ) then
593 
594  do j = 1, self%PureEFTGamma6%parameter_number
595  if ( i-n5 == j ) call self%PureEFTGamma6%parameter_names( j, name )
596  end do
597  return
598 
599  end if
600 
601  end subroutine eftcambpureeftstdparameternames
602 
603  ! ---------------------------------------------------------------------------------------------
605  subroutine eftcambpureeftstdparameternameslatex( self, i, latexname )
606 
607  implicit none
608 
609  class(eftcamb_std_pure_eft) :: self
610  integer , intent(in) :: i
611  character(*), intent(out) :: latexname
612 
613  integer :: nomega, nw, n1, n2, n3, n4, n5, n6
614  integer :: j
615 
616  ! compute the incremental number of parameters:
617  nomega = self%PureEFTOmega%parameter_number
618  nw = nomega + self%PureEFTwDE%parameter_number
619  n1 = nw + self%PureEFTGamma1%parameter_number
620  n2 = n1 + self%PureEFTGamma2%parameter_number
621  n3 = n2 + self%PureEFTGamma3%parameter_number
622  if ( .not. self%PureEFTHorndeski ) then
623  n4 = n3 + self%PureEFTGamma4%parameter_number
624  n5 = n4 + self%PureEFTGamma5%parameter_number
625  n6 = n5 + self%PureEFTGamma6%parameter_number
626  end if
627 
628  ! check validity of input:
629  if ( i > self%parameter_number .or. i <= 0) then
630  write(*,'(a,I3)') 'No parameter corresponding to: ', i
631  write(*,'(a,I3)') 'Total number of parameters is: ', self%parameter_number
632  call mpistop('EFTCAMB error')
633 
634  ! parameter from Omega function
635  else if ( i <= nomega ) then
636  do j = 1, self%PureEFTOmega%parameter_number
637  if ( i == j ) call self%PureEFTOmega%parameter_names_latex( j, latexname )
638  end do
639  return
640 
641  ! parameter from wDE function
642  else if ( i <= nw ) then
643  do j = 1, self%PureEFTwDE%parameter_number
644  if ( i-nomega == j ) call self%PureEFTwDE%parameter_names_latex( j, latexname )
645  end do
646  return
647 
648  ! parameter from Gamma1 function
649  else if ( i <= n1) then
650  do j = 1, self%PureEFTGamma1%parameter_number
651  if ( i-nw == j ) call self%PureEFTGamma1%parameter_names_latex( j, latexname )
652  end do
653  return
654 
655  !parameter from Gamma2 function
656  else if ( i <= n2 ) then
657  do j = 1, self%PureEFTGamma2%parameter_number
658  if ( i-n1 == j ) call self%PureEFTGamma2%parameter_names_latex( j, latexname )
659  end do
660  return
661 
662  !parameter from Gamma3 function
663  else if ( i <= n3 ) then
664  do j = 1, self%PureEFTGamma3%parameter_number
665  if ( i-n2 == j ) call self%PureEFTGamma3%parameter_names_latex( j, latexname )
666  end do
667  return
668 
669  !parameter from Gamma4 function
670  else if ( .not. self%PureEFTHorndeski .and. i <= n4 ) then
671  do j = 1, self%PureEFTGamma4%parameter_number
672  if ( i-n3 == j ) call self%PureEFTGamma4%parameter_names_latex( j, latexname )
673  end do
674  return
675 
676  !parameter from Gamma5 function
677  else if ( .not. self%PureEFTHorndeski .and. i <= n5 ) then
678  do j = 1, self%PureEFTGamma5%parameter_number
679  if ( i-n4 == j ) call self%PureEFTGamma5%parameter_names_latex( j, latexname )
680  end do
681  return
682 
683  !parameter from Gamma6 function
684  else if ( .not. self%PureEFTHorndeski .and. i <= n6 ) then
685  do j = 1, self%PureEFTGamma6%parameter_number
686  if ( i-n5 == j ) call self%PureEFTGamma6%parameter_names_latex( j, latexname )
687  end do
688  return
689 
690  end if
691 
692  end subroutine eftcambpureeftstdparameternameslatex
693 
694  ! ---------------------------------------------------------------------------------------------
696  subroutine eftcambpureeftstdparametervalues( self, i, value )
697 
698  implicit none
699 
700  class(eftcamb_std_pure_eft) :: self
701  integer , intent(in) :: i
702  real(dl), intent(out) :: value
703 
704  integer :: nomega, nw, n1, n2, n3, n4, n5, n6
705  integer :: j
706 
707  ! compute the incremental number of parameters:
708  nomega = self%PureEFTOmega%parameter_number
709  nw = nomega + self%PureEFTwDE%parameter_number
710  n1 = nw + self%PureEFTGamma1%parameter_number
711  n2 = n1 + self%PureEFTGamma2%parameter_number
712  n3 = n2 + self%PureEFTGamma3%parameter_number
713  if ( .not. self%PureEFTHorndeski ) then
714  n4 = n3 + self%PureEFTGamma4%parameter_number
715  n5 = n4 + self%PureEFTGamma5%parameter_number
716  n6 = n5 + self%PureEFTGamma6%parameter_number
717  end if
718 
719  ! check validity of input:
720  if ( i > self%parameter_number .or. i <= 0) then
721  write(*,'(a,I3)') 'EFTCAMB error: no parameter corresponding to number ',i
722  write(*,'(a,I3)') 'Total number of parameters is ', self%parameter_number
723  call mpistop('EFTCAMB error')
724 
725  ! parameter from Omega function
726  else if ( i <= nomega ) then
727  do j = 1, self%PureEFTOmega%parameter_number
728  if ( i == j ) call self%PureEFTOmega%parameter_value( j, value )
729  end do
730  return
731 
732  ! parameter from wDE function
733  else if ( i <= nw ) then
734  do j = 1, self%PureEFTwDE%parameter_number
735  if ( i-nomega == j ) call self%PureEFTwDE%parameter_value( j, value )
736  end do
737  return
738 
739  ! parameter from Gamma1 function
740  else if ( i <= n1) then
741  do j = 1, self%PureEFTGamma1%parameter_number
742  if ( i-nw == j ) call self%PureEFTGamma1%parameter_value( j, value )
743  end do
744  return
745 
746  !parameter from Gamma2 function
747  else if ( i <= n2 ) then
748  do j = 1, self%PureEFTGamma2%parameter_number
749  if ( i-n1 == j ) call self%PureEFTGamma2%parameter_value( j, value )
750  end do
751  return
752 
753  !parameter from Gamma3 function
754  else if ( i <= n3 ) then
755  do j = 1, self%PureEFTGamma3%parameter_number
756  if ( i-n2 == j ) call self%PureEFTGamma3%parameter_value( j, value )
757  end do
758  return
759 
760  !parameter from Gamma4 function
761  else if ( .not. self%PureEFTHorndeski .and. i <= n4 ) then
762  do j = 1, self%PureEFTGamma4%parameter_number
763  if ( i-n3 == j ) call self%PureEFTGamma4%parameter_value( j, value )
764  end do
765  return
766 
767  !parameter from Gamma5 function
768  else if ( .not. self%PureEFTHorndeski .and. i <= n5 ) then
769  do j = 1, self%PureEFTGamma5%parameter_number
770  if ( i-n4 == j ) call self%PureEFTGamma5%parameter_value( j, value )
771  end do
772  return
773 
774  !parameter from Gamma6 function
775  else if ( .not. self%PureEFTHorndeski .and. i <= n6 ) then
776  do j = 1, self%PureEFTGamma6%parameter_number
777  if ( i-n5 == j ) call self%PureEFTGamma6%parameter_value( j, value )
778  end do
779  return
780 
781  end if
782 
783  end subroutine eftcambpureeftstdparametervalues
784 
785  ! ---------------------------------------------------------------------------------------------
787  subroutine eftcambpureeftstdbackgroundeftfunctions( self, a, eft_par_cache, eft_cache )
788 
789  implicit none
790 
791  class(eftcamb_std_pure_eft) :: self
792  real(dl), intent(in) :: a
793  type(eftcamb_parameter_cache), intent(inout) :: eft_par_cache
794  type(eftcamb_timestep_cache ), intent(inout) :: eft_cache
795 
796  eft_cache%EFTOmegaV = self%PureEFTOmega%value(a)
797  eft_cache%EFTOmegaP = self%PureEFTOmega%first_derivative(a)
798  eft_cache%EFTOmegaPP = self%PureEFTOmega%second_derivative(a)
799  eft_cache%EFTOmegaPPP = self%PureEFTOmega%third_derivative(a)
800  eft_cache%EFTc = ( eft_cache%adotoa**2 - eft_cache%Hdot )*( eft_cache%EFTOmegaV + 0.5_dl*a*eft_cache%EFTOmegaP ) &
801  & -0.5_dl*( a*eft_cache%adotoa )**2*eft_cache%EFTOmegaPP&
802  & +0.5_dl*eft_cache%grhov_t*( 1._dl+ self%PureEFTwDE%value(a) )
803  eft_cache%EFTLambda = -eft_cache%EFTOmegaV*( 2._dl*eft_cache%Hdot +eft_cache%adotoa**2 ) &
804  & -a*eft_cache%EFTOmegaP*( 2._dl*eft_cache%adotoa**2 + eft_cache%Hdot ) &
805  & -( a*eft_cache%adotoa )**2*eft_cache%EFTOmegaPP &
806  & +self%PureEFTwDE%value(a)*eft_cache%grhov_t
807  eft_cache%EFTcdot = +0.5_dl*eft_cache%adotoa*eft_cache%grhov_t*( -3._dl*(1._dl +self%PureEFTwDE%value(a))**2 + a*self%PureEFTwDE%first_derivative(a) ) &
808  & -eft_cache%EFTOmegaV*( eft_cache%Hdotdot -4._dl*eft_cache%adotoa*eft_cache%Hdot +2._dl*eft_cache%adotoa**3 ) &
809  & +0.5_dl*a*eft_cache%EFTOmegaP*( -eft_cache%Hdotdot +eft_cache%adotoa*eft_cache%Hdot +eft_cache%adotoa**3) &
810  & +0.5_dl*a**2*eft_cache%adotoa*eft_cache%EFTOmegaPP*( eft_cache%adotoa**2 -3._dl*eft_cache%Hdot ) &
811  & -0.5_dl*(a*eft_cache%adotoa)**3*eft_cache%EFTOmegaPPP
812  eft_cache%EFTLambdadot = -2._dl*eft_cache%EFTOmegaV*( eft_cache%Hdotdot -eft_cache%adotoa*eft_cache%Hdot -eft_cache%adotoa**3 ) &
813  & -a*eft_cache%EFTOmegaP*( +eft_cache%Hdotdot +5._dl*eft_cache%adotoa*eft_cache%Hdot -eft_cache%adotoa**3 ) &
814  & -a**2*eft_cache%EFTOmegaPP*eft_cache%adotoa*( +2._dl*eft_cache%adotoa**2 +3._dl*eft_cache%Hdot )&
815  & -(a*eft_cache%adotoa)**3*eft_cache%EFTOmegaPPP &
816  & +eft_cache%grhov_t*eft_cache%adotoa*( a*self%PureEFTwDE%first_derivative(a) -3._dl*self%PureEFTwDE%value(a)*(1._dl +self%PureEFTwDE%value(a) ))
817 
818  end subroutine eftcambpureeftstdbackgroundeftfunctions
819 
820  ! ---------------------------------------------------------------------------------------------
822  subroutine eftcambpureeftstdsecondordereftfunctions( self, a, eft_par_cache, eft_cache )
823 
824  implicit none
825 
826  class(eftcamb_std_pure_eft) :: self
827  real(dl), intent(in) :: a
828  type(eftcamb_parameter_cache), intent(inout) :: eft_par_cache
829  type(eftcamb_timestep_cache ), intent(inout) :: eft_cache
830 
831  eft_cache%EFTGamma1V = self%PureEFTGamma1%value(a)
832  eft_cache%EFTGamma1P = self%PureEFTGamma1%first_derivative(a)
833  eft_cache%EFTGamma2V = self%PureEFTGamma2%value(a)
834  eft_cache%EFTGamma2P = self%PureEFTGamma2%first_derivative(a)
835  eft_cache%EFTGamma3V = self%PureEFTGamma3%value(a)
836  eft_cache%EFTGamma3P = self%PureEFTGamma3%first_derivative(a)
837  if ( self%PureEFTHorndeski ) then
838  eft_cache%EFTGamma4V = -eft_cache%EFTGamma3V
839  eft_cache%EFTGamma4P = -eft_cache%EFTGamma3P
840  eft_cache%EFTGamma4PP = -self%PureEFTGamma3%second_derivative(a)
841  eft_cache%EFTGamma5V = +0.5_dl*eft_cache%EFTGamma3V
842  eft_cache%EFTGamma5P = +0.5_dl*eft_cache%EFTGamma3P
843  eft_cache%EFTGamma6V = 0._dl
844  eft_cache%EFTGamma6P = 0._dl
845  else
846  eft_cache%EFTGamma4V = self%PureEFTGamma4%value(a)
847  eft_cache%EFTGamma4P = self%PureEFTGamma4%first_derivative(a)
848  eft_cache%EFTGamma4PP = self%PureEFTGamma4%second_derivative(a)
849  eft_cache%EFTGamma5V = self%PureEFTGamma5%value(a)
850  eft_cache%EFTGamma5P = self%PureEFTGamma5%first_derivative(a)
851  eft_cache%EFTGamma6V = self%PureEFTGamma6%value(a)
852  eft_cache%EFTGamma6P = self%PureEFTGamma6%first_derivative(a)
853  end if
854 
855  end subroutine eftcambpureeftstdsecondordereftfunctions
856 
857  ! ---------------------------------------------------------------------------------------------
860  function eftcambpureeftstdcomputedtauda( self, a, eft_par_cache, eft_cache )
861 
862  implicit none
863 
864  class(eftcamb_std_pure_eft) :: self
865  real(dl), intent(in) :: a
866  type(eftcamb_parameter_cache), intent(inout) :: eft_par_cache
867  type(eftcamb_timestep_cache ), intent(inout) :: eft_cache
868 
869  real(dl) :: eftcambpureeftstdcomputedtauda
870 
871  real(dl) :: temp
872 
873  temp = eft_cache%grhoa2 +eft_par_cache%grhov*a*a*self%PureEFTwDE%integral(a)
874  eftcambpureeftstdcomputedtauda = sqrt(3/temp)
875 
876  end function eftcambpureeftstdcomputedtauda
877 
878  ! ---------------------------------------------------------------------------------------------
880  subroutine eftcambpureeftstdcomputeadotoa( self, a, eft_par_cache, eft_cache )
881 
882  implicit none
883 
884  class(eftcamb_std_pure_eft) :: self
885  real(dl), intent(in) :: a
886  type(eftcamb_parameter_cache), intent(inout) :: eft_par_cache
887  type(eftcamb_timestep_cache ), intent(inout) :: eft_cache
888 
889  eft_cache%grhov_t = eft_par_cache%grhov*self%PureEFTwDE%integral(a)
890  eft_cache%adotoa = sqrt( ( eft_cache%grhom_t +eft_cache%grhov_t )/3._dl )
891 
892  end subroutine eftcambpureeftstdcomputeadotoa
893 
894  ! ---------------------------------------------------------------------------------------------
896  subroutine eftcambpureeftstdcomputehubbleder( self, a, eft_par_cache, eft_cache )
897 
898  implicit none
899 
900  class(eftcamb_std_pure_eft) :: self
901  real(dl), intent(in) :: a
902  type(eftcamb_parameter_cache), intent(inout) :: eft_par_cache
903  type(eftcamb_timestep_cache ), intent(inout) :: eft_cache
904 
905  eft_cache%gpiv_t = self%PureEFTwDE%value(a)*eft_cache%grhov_t
906  eft_cache%Hdot = -0.5_dl*( eft_cache%adotoa**2 +eft_cache%gpresm_t +eft_cache%gpiv_t )
907  eft_cache%Hdotdot = eft_cache%adotoa*( ( eft_cache%grhob_t +eft_cache%grhoc_t)/6._dl +2._dl*( eft_cache%grhor_t +eft_cache%grhog_t)/3._dl ) &
908  & +eft_cache%adotoa*eft_cache%grhov_t*( 1._dl/6._dl +self%PureEFTwDE%value(a) +1.5_dl*self%PureEFTwDE%value(a)**2 -0.5_dl*a*self%PureEFTwDE%first_derivative(a) ) &
909  & +eft_cache%adotoa*eft_cache%grhonu_tot/6._dl -0.5_dl*eft_cache%adotoa*eft_cache%gpinu_tot -0.5_dl*eft_cache%gpinudot_tot
910 
911  end subroutine eftcambpureeftstdcomputehubbleder
912 
913  ! ---------------------------------------------------------------------------------------------
915  function eftcambpureeftstdadditionalmodelstability( self, a, eft_par_cache, eft_cache )
916 
917  implicit none
918 
919  class(eftcamb_std_pure_eft) :: self
920  real(dl), intent(in) :: a
921  type(eftcamb_parameter_cache), intent(inout) :: eft_par_cache
922  type(eftcamb_timestep_cache ), intent(inout) :: eft_cache
923 
924  logical :: eftcambpureeftstdadditionalmodelstability
925 
926  eftcambpureeftstdadditionalmodelstability = .true.
927  if ( self%PureEFTwDE%value(a) > -1._dl/3._dl ) eftcambpureeftstdadditionalmodelstability = .false.
928 
929  end function eftcambpureeftstdadditionalmodelstability
930 
931  ! ---------------------------------------------------------------------------------------------
932 
933 end module eftcamb_pure_eft_std
934 
935 !----------------------------------------------------------------------------------------
This module contains the definition of the constant parametrization, inheriting from parametrized_fun...
This module contains the definition of the Taylor expansion parametrization, around a=0...
This module contains the definition of the EFTCAMB caches. These are used to store parameters that ca...
This module contains the definition of the linear parametrization, inheriting from parametrized_funct...
This module contains the definition of the turning point parametrization, inheriting from parametrize...
This module contains the abstract definition of all the places where EFTCAMB interacts with CAMB in c...
This module contains the definition of the CPL parametrization, inheriting from parametrized_function...
This module contains the definition of the power law parametrization, inheriting from parametrized_fu...
This module contains the definition of neutral parametrizations that can be used when parametrized fu...
This module contains the definitions of all the EFTCAMB compile time flags.
Definition: 01_EFT_def.f90:25
This module contains the definition of the generalized Jassal-Bagla-Padmanabhan (JBP) parametrization...
This module contains the definition of the Pure EFT model in which the EFT is described by six functi...
This module contains the definition of the exponential parametrization, inheriting from parametrized_...
This module contains the abstract class for generic parametrizations for 1D functions that are used b...