insolation¶
digraph inheritance1026093583 { bgcolor=transparent; rankdir=LR; ratio=expand; size=""; "AnnualMeanInsolation" [URL="#climlab.radiation.insolation.AnnualMeanInsolation",dirType=back,fillcolor=white,fontname="Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans",fontsize=14,height=0.25,shape=box,style="setlinewidth(0.5),filled",target="_top",tooltip="A class for latitudewise solar insolation averaged over a year."]; "_Insolation" -> "AnnualMeanInsolation" [arrowsize=0.5,dirType=back,style="setlinewidth(0.5)"]; "DailyInsolation" [URL="#climlab.radiation.insolation.DailyInsolation",dirType=back,fillcolor=white,fontname="Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans",fontsize=14,height=0.25,shape=box,style="setlinewidth(0.5),filled",target="_top",tooltip="A class to compute latitudewise daily solar insolation for specific"]; "AnnualMeanInsolation" -> "DailyInsolation" [arrowsize=0.5,dirType=back,style="setlinewidth(0.5)"]; "DiagnosticProcess" [URL="climlab.process.diagnostic.html#climlab.process.diagnostic.DiagnosticProcess",dirType=back,fillcolor=white,fontname="Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans",fontsize=14,height=0.25,shape=box,style="setlinewidth(0.5),filled",target="_top",tooltip="A parent class for all processes that are strictly diagnostic,"]; "TimeDependentProcess" -> "DiagnosticProcess" [arrowsize=0.5,dirType=back,style="setlinewidth(0.5)"]; "FixedInsolation" [URL="#climlab.radiation.insolation.FixedInsolation",dirType=back,fillcolor=white,fontname="Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans",fontsize=14,height=0.25,shape=box,style="setlinewidth(0.5),filled",target="_top",tooltip="A class for fixed insolation at each point of latitude off the domain."]; "_Insolation" -> "FixedInsolation" [arrowsize=0.5,dirType=back,style="setlinewidth(0.5)"]; "InstantInsolation" [URL="#climlab.radiation.insolation.InstantInsolation",dirType=back,fillcolor=white,fontname="Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans",fontsize=14,height=0.25,shape=box,style="setlinewidth(0.5),filled",target="_top",tooltip="A class to compute latitudewise instantaneous solar insolation for specific"]; "AnnualMeanInsolation" -> "InstantInsolation" [arrowsize=0.5,dirType=back,style="setlinewidth(0.5)"]; "P2Insolation" [URL="#climlab.radiation.insolation.P2Insolation",dirType=back,fillcolor=white,fontname="Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans",fontsize=14,height=0.25,shape=box,style="setlinewidth(0.5),filled",target="_top",tooltip="A class for parabolic solar distribution over the domain's latitude"]; "_Insolation" -> "P2Insolation" [arrowsize=0.5,dirType=back,style="setlinewidth(0.5)"]; "Process" [URL="climlab.process.process.html#climlab.process.process.Process",dirType=back,fillcolor=white,fontname="Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans",fontsize=14,height=0.25,shape=box,style="setlinewidth(0.5),filled",target="_top",tooltip="A generic parent class for all climlab process objects."]; "TimeDependentProcess" [URL="climlab.process.time_dependent_process.html#climlab.process.time_dependent_process.TimeDependentProcess",dirType=back,fillcolor=white,fontname="Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans",fontsize=14,height=0.25,shape=box,style="setlinewidth(0.5),filled",target="_top",tooltip="A generic parent class for all time-dependent processes."]; "Process" -> "TimeDependentProcess" [arrowsize=0.5,dirType=back,style="setlinewidth(0.5)"]; "_Insolation" [dirType=back,fillcolor=white,fontname="Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans",fontsize=14,height=0.25,shape=box,style="setlinewidth(0.5),filled",tooltip="A private parent class for insolation processes."]; "DiagnosticProcess" -> "_Insolation" [arrowsize=0.5,dirType=back,style="setlinewidth(0.5)"]; }Classes to provide insolation as input for other CLIMLAB processes.
Options include
climlab.radiation.P2Insolation
(idealized 2nd Legendre polynomial form)climlab.radiation.FixedInsolation
(generic steady-state insolation)climlab.radiation.AnnualMeanInsolation
(steady-state annual-mean insolation computed from orbital parameters and latitude)climlab.radiation.DailyInsolation
(time-varying daily-mean insolation computed from orbital parameters, latitude and time of year)
All are subclasses of climlab.process.DiagnosticProcess
and do add any tendencies to any state variables.
At least two diagnostics are provided:
insolation
, the incoming solar radiation in \(\textrm{W}}{\textrm{m}^2}`\)coszen
, cosine of the solar zenith angle
- class climlab.radiation.insolation.AnnualMeanInsolation(S0=1365.2, orb={'ecc': 0.017236, 'long_peri': 281.37, 'obliquity': 23.446}, **kwargs)[source]¶
Bases:
_Insolation
A class for latitudewise solar insolation averaged over a year.
This class computes the solar insolation for each day of the year and latitude specified in the domain on the basis of orbital parameters and astronomical formulas.
Therefore it uses the method
daily_insolation()
. For details how the solar distribution is dependend on orbital parameters see there.The mean over the year is calculated from data given by
daily_insolation()
and stored in the object’s attributeself.insolation
Initialization parameters
- Parameters:
S0 (float) –
solar constant
unit: \(\frac{\textrm{W}}{\textrm{m}^2}\)
default value:
1365.2
orb (dict) –
a dictionary with three orbital parameters (as provided by
OrbitalTable
):'ecc'
- eccentricityunit: dimensionless
default value:
0.017236
'long_peri'
- longitude of perihelion (precession angle)unit: degrees
default value:
281.37
'obliquity'
- obliquity angleunit: degrees
default value:
23.446
Object attributes
Additional to the parent class
_Insolation
following object attributes are generated and updated during initialization:- Variables:
- Example:
Create regular EBM and replace standard insolation subprocess by
AnnualMeanInsolation
:>>> import climlab >>> from climlab.radiation import AnnualMeanInsolation >>> # model creation >>> model = climlab.EBM() >>> print model
climlab Process of type <class 'climlab.model.ebm.EBM'>. State variables and domain shapes: Ts: (90, 1) The subprocess tree: top: <class 'climlab.model.ebm.EBM'> diffusion: <class 'climlab.dynamics.diffusion.MeridionalDiffusion'> LW: <class 'climlab.radiation.AplusBT.AplusBT'> albedo: <class 'climlab.surface.albedo.StepFunctionAlbedo'> iceline: <class 'climlab.surface.albedo.Iceline'> cold_albedo: <class 'climlab.surface.albedo.ConstantAlbedo'> warm_albedo: <class 'climlab.surface.albedo.P2Albedo'> insolation: <class 'climlab.radiation.insolation.P2Insolation'>
>>> # catch model domain for subprocess creation >>> sfc = model.domains['Ts'] >>> # create AnnualMeanInsolation subprocess >>> new_insol = AnnualMeanInsolation(domains=sfc, **model.param) >>> # add it to the model >>> model.add_subprocess('insolation',new_insol) >>> print model
climlab Process of type <class 'climlab.model.ebm.EBM'>. State variables and domain shapes: Ts: (90, 1) The subprocess tree: top: <class 'climlab.model.ebm.EBM'> diffusion: <class 'climlab.dynamics.diffusion.MeridionalDiffusion'> LW: <class 'climlab.radiation.AplusBT.AplusBT'> albedo: <class 'climlab.surface.albedo.StepFunctionAlbedo'> iceline: <class 'climlab.surface.albedo.Iceline'> cold_albedo: <class 'climlab.surface.albedo.ConstantAlbedo'> warm_albedo: <class 'climlab.surface.albedo.P2Albedo'> insolation: <class 'climlab.radiation.insolation.AnnualMeanInsolation'>
- Attributes:
S0
Property of solar constant S0.
depth
Depth at grid centers (m)
depth_bounds
Depth at grid interfaces (m)
diagnostics
Dictionary access to all diagnostic variables
input
Dictionary access to all input variables
lat
Latitude of grid centers (degrees North)
lat_bounds
Latitude of grid interfaces (degrees North)
lev
Pressure levels at grid centers (hPa or mb)
lev_bounds
Pressure levels at grid interfaces (hPa or mb)
lon
Longitude of grid centers (degrees)
lon_bounds
Longitude of grid interfaces (degrees)
orb
Property of dictionary for orbital parameters.
timestep
The amount of time over which
step_forward()
is integrating in unit seconds.
Methods
add_diagnostic
(name[, value])Create a new diagnostic variable called
name
for this process and initialize it with the givenvalue
.add_input
(name[, value])Create a new input variable called
name
for this process and initialize it with the givenvalue
.add_subprocess
(name, proc)Adds a single subprocess to this process.
add_subprocesses
(procdict)Adds a dictionary of subproceses to this process.
compute
()Computes the tendencies for all state variables given current state and specified input.
compute_diagnostics
([num_iter])Compute all tendencies and diagnostics, but don't update model state.
declare_diagnostics
(diaglist)Add the variable names in
inputlist
to the list of diagnostics.declare_input
(inputlist)Add the variable names in
inputlist
to the list of necessary inputs.integrate_converge
([crit, verbose])Integrates the model until model states are converging.
integrate_days
([days, verbose])Integrates the model forward for a specified number of days.
integrate_years
([years, verbose])Integrates the model by a given number of years.
remove_diagnostic
(name)Removes a diagnostic from the
process.diagnostic
dictionary and also delete the associated process attribute.remove_subprocess
(name[, verbose])Removes a single subprocess from this process.
set_state
(name, value)Sets the variable
name
to a new statevalue
.set_timestep
([timestep, num_steps_per_year])Calculates the timestep in unit seconds and calls the setter function of
timestep()
step_forward
()Updates state variables with computed tendencies.
to_xarray
([diagnostics])Convert process variables to
xarray.Dataset
format.- property orb¶
Property of dictionary for orbital parameters.
orb contains: (for more information see
OrbitalTable
)'ecc'
- eccentricity [unit: dimensionless]'long_peri'
- longitude of perihelion (precession angle) [unit: degrees]'obliquity'
- obliquity angle [unit: degrees]
- Getter:
Returns the orbital dictionary which is stored in attribute
self._orb
.- Setter:
sets orb which is addressed as
self._orb
to the new valueupdates the parameter dictionary
self.param['orb']
andcalls method
_compute_fixed()
- Type:
- class climlab.radiation.insolation.DailyInsolation(S0=1365.2, orb={'ecc': 0.017236, 'long_peri': 281.37, 'obliquity': 23.446}, **kwargs)[source]¶
Bases:
AnnualMeanInsolation
A class to compute latitudewise daily solar insolation for specific days of the year.
This class computes the solar insolation on basis of orbital parameters and astronomical formulas.
Therefore it uses the method
daily_insolation()
. For details how the solar distribution is dependend on orbital parameters see there.Initialization parameters
- Parameters:
S0 (float) –
solar constant
unit: \(\frac{\textrm{W}}{\textrm{m}^2}\)
default value:
1365.2
orb (dict) –
a dictionary with orbital parameters:
'ecc'
- eccentricityunit: dimensionless
default value:
0.017236
'long_peri'
- longitude of perihelion (precession angle)unit: degrees
default value:
281.37
'obliquity'
- obliquity angleunit: degrees
default value:
23.446
Object attributes
Additional to the parent class
_Insolation
following object attributes are generated and updated during initialization:- Variables:
- Example:
Create regular EBM and replace standard insolation subprocess by
DailyInsolation
:>>> import climlab >>> from climlab.radiation import DailyInsolation >>> # model creation >>> model = climlab.EBM() >>> print model
climlab Process of type <class 'climlab.model.ebm.EBM'>. State variables and domain shapes: Ts: (90, 1) The subprocess tree: top: <class 'climlab.model.ebm.EBM'> diffusion: <class 'climlab.dynamics.diffusion.MeridionalDiffusion'> LW: <class 'climlab.radiation.AplusBT.AplusBT'> albedo: <class 'climlab.surface.albedo.StepFunctionAlbedo'> iceline: <class 'climlab.surface.albedo.Iceline'> cold_albedo: <class 'climlab.surface.albedo.ConstantAlbedo'> warm_albedo: <class 'climlab.surface.albedo.P2Albedo'> insolation: <class 'climlab.radiation.insolation.P2Insolation'>
>>> # catch model domain for subprocess creation >>> sfc = model.domains['Ts'] >>> # create DailyInsolation subprocess and add it to the model >>> model.add_subprocess('insolation',DailyInsolation(domains=sfc, **model.param)) >>> print model
climlab Process of type <class 'climlab.model.ebm.EBM'>. State variables and domain shapes: Ts: (90, 1) The subprocess tree: top: <class 'climlab.model.ebm.EBM'> diffusion: <class 'climlab.dynamics.diffusion.MeridionalDiffusion'> LW: <class 'climlab.radiation.AplusBT.AplusBT'> albedo: <class 'climlab.surface.albedo.StepFunctionAlbedo'> iceline: <class 'climlab.surface.albedo.Iceline'> cold_albedo: <class 'climlab.surface.albedo.ConstantAlbedo'> warm_albedo: <class 'climlab.surface.albedo.P2Albedo'> insolation: <class 'climlab.radiation.insolation.DailyInsolation'>
- Attributes:
S0
Property of solar constant S0.
depth
Depth at grid centers (m)
depth_bounds
Depth at grid interfaces (m)
diagnostics
Dictionary access to all diagnostic variables
input
Dictionary access to all input variables
lat
Latitude of grid centers (degrees North)
lat_bounds
Latitude of grid interfaces (degrees North)
lev
Pressure levels at grid centers (hPa or mb)
lev_bounds
Pressure levels at grid interfaces (hPa or mb)
lon
Longitude of grid centers (degrees)
lon_bounds
Longitude of grid interfaces (degrees)
orb
Property of dictionary for orbital parameters.
timestep
The amount of time over which
step_forward()
is integrating in unit seconds.
Methods
add_diagnostic
(name[, value])Create a new diagnostic variable called
name
for this process and initialize it with the givenvalue
.add_input
(name[, value])Create a new input variable called
name
for this process and initialize it with the givenvalue
.add_subprocess
(name, proc)Adds a single subprocess to this process.
add_subprocesses
(procdict)Adds a dictionary of subproceses to this process.
compute
()Computes the tendencies for all state variables given current state and specified input.
compute_diagnostics
([num_iter])Compute all tendencies and diagnostics, but don't update model state.
declare_diagnostics
(diaglist)Add the variable names in
inputlist
to the list of diagnostics.declare_input
(inputlist)Add the variable names in
inputlist
to the list of necessary inputs.integrate_converge
([crit, verbose])Integrates the model until model states are converging.
integrate_days
([days, verbose])Integrates the model forward for a specified number of days.
integrate_years
([years, verbose])Integrates the model by a given number of years.
remove_diagnostic
(name)Removes a diagnostic from the
process.diagnostic
dictionary and also delete the associated process attribute.remove_subprocess
(name[, verbose])Removes a single subprocess from this process.
set_state
(name, value)Sets the variable
name
to a new statevalue
.set_timestep
([timestep, num_steps_per_year])Calculates the timestep in unit seconds and calls the setter function of
timestep()
step_forward
()Updates state variables with computed tendencies.
to_xarray
([diagnostics])Convert process variables to
xarray.Dataset
format.
- class climlab.radiation.insolation.FixedInsolation(S0=341.3, **kwargs)[source]¶
Bases:
_Insolation
A class for fixed insolation at each point of latitude off the domain.
The solar distribution for the whole domain is constant and specified by a parameter.
Initialization parameters
- Parameters:
S0 (float) –
solar constant
unit: \(\frac{\textrm{W}}{\textrm{m}^2}\)
default value:
const.S0/4 = 341.2
- Example:
>>> import climlab >>> from climlab.radiation.insolation import FixedInsolation >>> model = climlab.EBM() >>> sfc = model.Ts.domain >>> fixed_ins = FixedInsolation(S0=340.0, domains=sfc) >>> print fixed_ins climlab Process of type <class 'climlab.radiation.insolation.FixedInsolation'>. State variables and domain shapes: The subprocess tree: top: <class 'climlab.radiation.insolation.FixedInsolation'>
- Attributes:
S0
Property of solar constant S0.
depth
Depth at grid centers (m)
depth_bounds
Depth at grid interfaces (m)
diagnostics
Dictionary access to all diagnostic variables
input
Dictionary access to all input variables
lat
Latitude of grid centers (degrees North)
lat_bounds
Latitude of grid interfaces (degrees North)
lev
Pressure levels at grid centers (hPa or mb)
lev_bounds
Pressure levels at grid interfaces (hPa or mb)
lon
Longitude of grid centers (degrees)
lon_bounds
Longitude of grid interfaces (degrees)
timestep
The amount of time over which
step_forward()
is integrating in unit seconds.
Methods
add_diagnostic
(name[, value])Create a new diagnostic variable called
name
for this process and initialize it with the givenvalue
.add_input
(name[, value])Create a new input variable called
name
for this process and initialize it with the givenvalue
.add_subprocess
(name, proc)Adds a single subprocess to this process.
add_subprocesses
(procdict)Adds a dictionary of subproceses to this process.
compute
()Computes the tendencies for all state variables given current state and specified input.
compute_diagnostics
([num_iter])Compute all tendencies and diagnostics, but don't update model state.
declare_diagnostics
(diaglist)Add the variable names in
inputlist
to the list of diagnostics.declare_input
(inputlist)Add the variable names in
inputlist
to the list of necessary inputs.integrate_converge
([crit, verbose])Integrates the model until model states are converging.
integrate_days
([days, verbose])Integrates the model forward for a specified number of days.
integrate_years
([years, verbose])Integrates the model by a given number of years.
remove_diagnostic
(name)Removes a diagnostic from the
process.diagnostic
dictionary and also delete the associated process attribute.remove_subprocess
(name[, verbose])Removes a single subprocess from this process.
set_state
(name, value)Sets the variable
name
to a new statevalue
.set_timestep
([timestep, num_steps_per_year])Calculates the timestep in unit seconds and calls the setter function of
timestep()
step_forward
()Updates state variables with computed tendencies.
to_xarray
([diagnostics])Convert process variables to
xarray.Dataset
format.
- class climlab.radiation.insolation.InstantInsolation(S0=1365.2, orb={'ecc': 0.017236, 'long_peri': 281.37, 'obliquity': 23.446}, **kwargs)[source]¶
Bases:
AnnualMeanInsolation
A class to compute latitudewise instantaneous solar insolation for specific days of the year.
This class computes the solar insolation on basis of orbital parameters and astronomical formulas.
Therefore it uses the method
instant_insolation()
. For details how the solar distribution is dependend on orbital parameters see there.Initialization parameters
- Parameters:
S0 (float) –
solar constant
unit: \(\frac{\textrm{W}}{\textrm{m}^2}\)
default value:
1365.2
orb (dict) –
a dictionary with orbital parameters:
'ecc'
- eccentricityunit: dimensionless
default value:
0.017236
'long_peri'
- longitude of perihelion (precession angle)unit: degrees
default value:
281.37
'obliquity'
- obliquity angleunit: degrees
default value:
23.446
Object attributes
Additional to the parent class
_Insolation
following object attributes are generated and updated during initialization:- Variables:
- Example:
Create regular EBM and replace standard insolation subprocess by
DailyInsolation
:>>> import climlab >>> from climlab.radiation import InstantInsolation >>> # model creation >>> model = climlab.EBM() >>> print model
climlab Process of type <class 'climlab.model.ebm.EBM'>. State variables and domain shapes: Ts: (90, 1) The subprocess tree: top: <class 'climlab.model.ebm.EBM'> diffusion: <class 'climlab.dynamics.diffusion.MeridionalDiffusion'> LW: <class 'climlab.radiation.AplusBT.AplusBT'> albedo: <class 'climlab.surface.albedo.StepFunctionAlbedo'> iceline: <class 'climlab.surface.albedo.Iceline'> cold_albedo: <class 'climlab.surface.albedo.ConstantAlbedo'> warm_albedo: <class 'climlab.surface.albedo.P2Albedo'> insolation: <class 'climlab.radiation.insolation.P2Insolation'>
>>> # catch model domain for subprocess creation >>> sfc = model.domains['Ts'] >>> # create InstantInsolation subprocess and add it to the model >>> model.add_subprocess('insolation',InstantInsolation(domains=sfc, **model.param)) >>> print model
climlab Process of type <class 'climlab.model.ebm.EBM'>. State variables and domain shapes: Ts: (90, 1) The subprocess tree: top: <class 'climlab.model.ebm.EBM'> diffusion: <class 'climlab.dynamics.diffusion.MeridionalDiffusion'> LW: <class 'climlab.radiation.AplusBT.AplusBT'> albedo: <class 'climlab.surface.albedo.StepFunctionAlbedo'> iceline: <class 'climlab.surface.albedo.Iceline'> cold_albedo: <class 'climlab.surface.albedo.ConstantAlbedo'> warm_albedo: <class 'climlab.surface.albedo.P2Albedo'> insolation: <class 'climlab.radiation.insolation.InstantInsolation'>
- Attributes:
S0
Property of solar constant S0.
depth
Depth at grid centers (m)
depth_bounds
Depth at grid interfaces (m)
diagnostics
Dictionary access to all diagnostic variables
input
Dictionary access to all input variables
lat
Latitude of grid centers (degrees North)
lat_bounds
Latitude of grid interfaces (degrees North)
lev
Pressure levels at grid centers (hPa or mb)
lev_bounds
Pressure levels at grid interfaces (hPa or mb)
lon
Longitude of grid centers (degrees)
lon_bounds
Longitude of grid interfaces (degrees)
orb
Property of dictionary for orbital parameters.
timestep
The amount of time over which
step_forward()
is integrating in unit seconds.
Methods
add_diagnostic
(name[, value])Create a new diagnostic variable called
name
for this process and initialize it with the givenvalue
.add_input
(name[, value])Create a new input variable called
name
for this process and initialize it with the givenvalue
.add_subprocess
(name, proc)Adds a single subprocess to this process.
add_subprocesses
(procdict)Adds a dictionary of subproceses to this process.
compute
()Computes the tendencies for all state variables given current state and specified input.
compute_diagnostics
([num_iter])Compute all tendencies and diagnostics, but don't update model state.
declare_diagnostics
(diaglist)Add the variable names in
inputlist
to the list of diagnostics.declare_input
(inputlist)Add the variable names in
inputlist
to the list of necessary inputs.integrate_converge
([crit, verbose])Integrates the model until model states are converging.
integrate_days
([days, verbose])Integrates the model forward for a specified number of days.
integrate_years
([years, verbose])Integrates the model by a given number of years.
remove_diagnostic
(name)Removes a diagnostic from the
process.diagnostic
dictionary and also delete the associated process attribute.remove_subprocess
(name[, verbose])Removes a single subprocess from this process.
set_state
(name, value)Sets the variable
name
to a new statevalue
.set_timestep
([timestep, num_steps_per_year])Calculates the timestep in unit seconds and calls the setter function of
timestep()
step_forward
()Updates state variables with computed tendencies.
to_xarray
([diagnostics])Convert process variables to
xarray.Dataset
format.
- class climlab.radiation.insolation.P2Insolation(S0=1365.2, s2=-0.48, **kwargs)[source]¶
Bases:
_Insolation
A class for parabolic solar distribution over the domain’s latitude on the basis of the second order Legendre Polynomial.
Calculates the latitude dependent solar distribution as
\[S(\varphi) = \frac{S_0}{4} \left( 1 + s_2 P_2(x) \right)\]where \(P_2(x) = \frac{1}{2} (3x^2 - 1)\) is the second order Legendre Polynomial and \(x=sin(\varphi)\).
Initialization parameters
- Parameters:
S0 (float) –
solar constant
unit: \(\frac{\textrm{W}}{\textrm{m}^2}\)
default value:
1365.2
s2 (floar) –
factor for second legendre polynominal term
default value:
-0.48
- Example:
>>> import climlab >>> from climlab.radiation.insolation import P2Insolation >>> model = climlab.EBM() >>> sfc = model.Ts.domain >>> p2_ins = P2Insolation(S0=340.0, s2=-0.5, domains=sfc) >>> print p2_ins climlab Process of type <class 'climlab.radiation.insolation.P2Insolation'>. State variables and domain shapes: The subprocess tree: top: <class 'climlab.radiation.insolation.P2Insolation'>
- Attributes:
S0
Property of solar constant S0.
depth
Depth at grid centers (m)
depth_bounds
Depth at grid interfaces (m)
diagnostics
Dictionary access to all diagnostic variables
input
Dictionary access to all input variables
lat
Latitude of grid centers (degrees North)
lat_bounds
Latitude of grid interfaces (degrees North)
lev
Pressure levels at grid centers (hPa or mb)
lev_bounds
Pressure levels at grid interfaces (hPa or mb)
lon
Longitude of grid centers (degrees)
lon_bounds
Longitude of grid interfaces (degrees)
s2
Property of second legendre polynomial factor s2.
timestep
The amount of time over which
step_forward()
is integrating in unit seconds.
Methods
add_diagnostic
(name[, value])Create a new diagnostic variable called
name
for this process and initialize it with the givenvalue
.add_input
(name[, value])Create a new input variable called
name
for this process and initialize it with the givenvalue
.add_subprocess
(name, proc)Adds a single subprocess to this process.
add_subprocesses
(procdict)Adds a dictionary of subproceses to this process.
compute
()Computes the tendencies for all state variables given current state and specified input.
compute_diagnostics
([num_iter])Compute all tendencies and diagnostics, but don't update model state.
declare_diagnostics
(diaglist)Add the variable names in
inputlist
to the list of diagnostics.declare_input
(inputlist)Add the variable names in
inputlist
to the list of necessary inputs.integrate_converge
([crit, verbose])Integrates the model until model states are converging.
integrate_days
([days, verbose])Integrates the model forward for a specified number of days.
integrate_years
([years, verbose])Integrates the model by a given number of years.
remove_diagnostic
(name)Removes a diagnostic from the
process.diagnostic
dictionary and also delete the associated process attribute.remove_subprocess
(name[, verbose])Removes a single subprocess from this process.
set_state
(name, value)Sets the variable
name
to a new statevalue
.set_timestep
([timestep, num_steps_per_year])Calculates the timestep in unit seconds and calls the setter function of
timestep()
step_forward
()Updates state variables with computed tendencies.
to_xarray
([diagnostics])Convert process variables to
xarray.Dataset
format.- property s2¶
Property of second legendre polynomial factor s2.
s2 in following equation:
\[S(\varphi) = \frac{S_0}{4} \left( 1 + s_2 P_2(x) \right)\]- Getter:
Returns the s2 parameter which is stored in attribute
self._s2
.- Setter:
sets s2 which is addressed as
self._S0
to the new valueupdates the parameter dictionary
self.param['s2']
andcalls method
_compute_fixed()
- Type: