fronts.solve_from_guess
- fronts.solve_from_guess(D, i, b, o_guess, guess, radial=False, max_nodes=1000, verbose=0)
Alternative solver for problems with a Dirichlet boundary condition.
Given a positive function D, scalars \(\theta_i\), \(\theta_b\) and \(o_b\), and coordinate unit vector \(\mathbf{\hat{r}}\), finds a function \(\theta\) of r and t such that:
\[\begin{split}\begin{cases} \dfrac{\partial\theta}{\partial t} = \nabla\cdot\left[D(\theta))\dfrac{\partial\theta}{\partial r} \mathbf{\hat{r}}\right ] & r>r_b(t),t>0\\ \theta(r, 0) = \theta_i & r>0 \\ \theta(r_b(t), t) = \theta_b & t>0 \\ r_b(t) = o_b\sqrt t \end{cases}\end{split}\]Alternative to the main
solve()
function. This function requires a starting mesh and guess of the solution. It is significantly less robust thansolve()
, and will fail to converge in many cases that the latter can easily handle (whether it converges will usually depend heavily on the problem, the starting mesh and the guess of the solution; it will raise a RuntimeError on failure). However, when it converges it is usually faster thansolve()
, which may be an advantage for some use cases. You should nonetheless prefersolve()
unless you have a particular use case for which you have found this function to be better.Possible use cases include refining a solution (note that
solve()
can do that too), optimization runs in which known solutions make good first approximations of solutions with similar parameters and every second of computing time counts, and in the implementation of other solving algorithms. In all these cases,solve()
should probably be used as a fallback for when this function fails.- Parameters:
D (callable or sympy.Expression or str or float) –
Callable that evaluates \(D\) and its derivatives, obtained from the
fronts.D
module or defined in the same manner—i.e.:D(theta)
evaluates and returns \(D\) attheta
D(theta, 1)
returns both the value of \(D\) and its first derivative attheta
D(theta, 2)
returns the value of \(D\), its first derivative, and its second derivative attheta
where
theta
may be a single float or a NumPy array.Alternatively, instead of a callable, the argument can be the expression of \(D\) in the form of a string or
sympy.Expression
with a single variable. In this case, the solver will differentiate and evaluate the expression as necessary.i (float) – Initial condition, \(\theta_i\).
b (float) – Imposed boundary value, \(\theta_b\).
o_guess (numpy.array_like, shape (n_guess,)) –
Starting mesh in terms of the Boltzmann variable. Must be strictly increasing.
o_guess[0]
is taken as the value of the parameter \(o_b\), which determines the behavior of the boundary. If zero, it implies that the boundary always exists at \(r=0\). It must be strictly positive if radial is not False. A non-zero value implies a moving boundary.On the other end,
o_guess[-1]
must be large enough to contain the solution to the semi-infinite problem.guess (float or numpy.array_like, shape (n_guess,)) – Starting guess of the solution at the points in o_guess. If a single value, the guess is assumed uniform.
radial ({False, 'cylindrical', 'polar', 'spherical'}, optional) –
Choice of coordinate unit vector \(\mathbf{\hat{r}}\). Must be one of the following:
False (default): \(\mathbf{\hat{r}}\) is any coordinate unit vector in rectangular (Cartesian) coordinates, or an axial unit vector in a cylindrical coordinate system
'cylindrical'
or'polar'
: \(\mathbf{\hat{r}}\) is the radial unit vector in a cylindrical or polar coordinate system'spherical'
: \(\mathbf{\hat{r}}\) is the radial unit vector in a spherical coordinate system
max_nodes (int, optional) – Maximum allowed number of mesh nodes.
verbose ({0, 1, 2}, optional) –
Level of algorithm’s verbosity. Must be one of the following:
0 (default): work silently
1: display a termination report
2: also display progress during iterations
- Returns:
solution – See
Solution
for a description of the solution object. Additional fields specific to this solver are included in the object:o (numpy.ndarray, shape (n,)) – Final solver mesh, in terms of the Boltzmann variable.
niter (int) – Number of iterations required to find the solution.
rms_residuals (numpy.ndarray, shape (n-1,)) – RMS values of the relative residuals over each mesh interval.
- Return type:
See also
Notes
This function works by transforming the partial differential equation with the Boltzmann transformation using
ode()
and then solving the resulting ODE with SciPy’s collocation-based boundary value problem solverscipy.integrate.solve_bvp()
and a two-point Dirichlet condition that matches the boundary and initial conditions of the problem. Upon that solver’s convergence, it runs a final check on whether the candidate solution also satisfies the semi-infinite condition (which implies \(d\theta/do\to0\) as \(o\to\infty\)).