<content> dsol </content>

3ff9934e-0071-4049-99d1-c9dc1e9d79f9 311 Alvaro Diaz Falconi 4 5 false false 0 0 decimal

Solving Differential equations

<content> dsol </content> <content> Utils </content>

M Cols C 0 : k 1 M cols range : 12 mat C 1 k el M k col : 21 line :

a b n pR a b a - n / 0 n range * + eval 11 line :

A D C B TriDiagSys n B length : v n 1 matrix : y v : 13 mat w A 1 el 1 -^: y 1 el B 1 el w * : 12 mat k 2 n range v k 1 - el C k 1 - el w * eval : w A k el D k el v k 1 - el * - (1 -^eval : y k el B k el D k el y k 1 - el * - (w * eval : 31 line for k n 1 - (1 range y k el y k el v k el y k 1 + el * - eval : for y 51 line :

x y pStem XY 03 matrix : S 02 matrix : 12 mat k 1 x length range XY XY x k el y k el o 13 mat stack : S S x k el 0 x k el y k el x k el 032 mat stack : 21 line for XY S 21 sys eval 31 line :

e# a# b# Equrep equrep(e# num2str, a# num2str, b# num2str) concat str2num 11 line :

x# a# Find n# 1 x# length range : f# x# n# el n# augment eval a# 1 findrows : (0 ≡ 0 f# 2 col eval if 21 line :

fs# x# feval fs# num2str (x# num2str) concat str2num 11 line :

fs# x# y# feval fs# num2str (x# num2str, y# num2str) concat str2num 11 line :

fs# x# y# z# feval fs# num2str (x# num2str, y# num2str, z# num2str) concat str2num 11 line :

<content> Options </content>

O# o# Options Ω# 01 matrix eval : n# 0 : 12 mat ω# o# sys2mat.1 num2str vectorize eval ω# ω# " strrep : n# n# 1 + : 12 mat g# ω# ≡ strsplit : length 1 > f# 1 O# rows range O# augment eval g# 1 el 2 findrows : (0 ≠ (O# f# 1 el 2 el g# 2 el : O# O# g# 1 el g# 2 el augment 11 O# cols range el stack eval : if Ω# Ω# g# 2 el stack eval : 21 line f# O# ω# 1 findrows : (0 ≠ (Ω# Ω# f# 2 el stack eval : Ω# Ω#'NaN stack eval : if if 21 line for Ω# rows 2 > o# cols 1 ≡ Ω# Ω# transpose if Ω# 1 el if 31 line :

<content> Example Options </content>

dsol solver Rkadapt steps 100 stepsmax 0 eps 10^{-7} 42 mat :

dsol dummy Options'NaN

dsol steps dummy eps 31 mat Options 100'NaN 10^{-7} 31 mat

dsol steps dummy eps 13 mat Options 100'NaN 10^{-7} 13 mat

dsol steps solver rkfixed ≡ eps 31 mat Options 100 rkfixed 10^{-7} 31 mat

dsol solver rkfixed steps 100 stepsmax 0 eps 10^{-7} 42 mat

dsol Jacobian J# ≡ Options J#

dsol solver rkfixed steps 100 stepsmax 0 eps 10^{-7} Jacobian J# 52 mat

Check inside a procedure

x f dsol dummy x ≡ Options dsol 21 line :

1 f solver rkfixed steps 100 stepsmax 0 eps 10^{-7} Jacobian J# dummy 162 mat

dsol solver rkfixed steps 100 stepsmax 0 eps 10^{-7} Jacobian J# 52 mat

lbvp.2(f(x),[a,b],M,n) solves the linear BVP y'' + p(x)y' + q(x)y = φ(x) in a < x < b with    f(x) = [p(x), q(x), φ(x)],  M = [α1 β1 c1; α2 β2 c2]  and     α1*y(a) + β1*y'(b) = c1   and  α2*y(a) + β2*y'(b) = c2using central differences. 1 f x M n lbvp.2 h x 2 el x 1 el - n / : X x 1 el 0 n range h * + : ν n 1 + : P 0 : Q 0 : Φ 0 : α1 β1 c1 α2 β2 c2 23 mat M : 17 mat r 1 ν range P r el Q r el Φ r el 13 mat X r el f eval : for A D C B 14 mat 2 - h 2^Q * + (2 h P * - 2 / (2 h P * + 2 / (h 2^Φ * (14 mat eval : A 1 el D 1 el C 1 el B 1 el 14 mat β1 0 ≡ 100 c1 α1 / 14 mat A 1 el 2 h * α1 β1 / * D 1 el * + (02 B 1 el 2 h * c1 β1 / * D 1 el * + (14 mat if eval : A ν el D ν el C ν el B ν el 14 mat β2 0 ≡ 100 c2 α2 / 14 mat A ν el 2 h * α2 β2 / * C ν el * - (20 B ν el 2 h * c2 β2 / * C ν el * - (14 mat if eval : Y A D C B TriDiagSys : Y' Y : j 2 n range : 13 mat Y' j el Y j 1 + el Y j 1 - el - 2 h * / : Y' 1 el β1 0 ≡ Y 2 el Y 1 el - h / c1 β1 / α1 β1 / Y 1 el * - if : Y' ν el β2 0 ≡ Y ν el Y n el - h / c2 β2 / α2 β2 / Y ν el * - if : 13 mat X Y Y' augment eval 81 line :

bvp.2(f(x),[a,b],M,n,ε) solves the non-linear BVP y'' = f(x,y,y') in a < x < b with M = [α1 β1 c1; α2 β2 c2] where     α1*y(a) + β1*y'(b) = c1   and  α2*y(a) + β2*y'(b) = c2using central differences. fs# x Y M n ε bvp.2 x# y# y'# f fs# x# y# y'# feval : h x 2 el x 1 el - n / : X x 1 el 0 n range h * + : ν n 1 + : P 0 : Q 0 : Φ 0 : α1 β1 c1 α2 β2 c2 23 mat M : 17 mat Y 0 ≡ ψ x 1 el α1 * β1 + α1 x 2 el α2 * β2 + α2 22 mat (1 -^c1 c1 21 mat * : Y ψ 2 el ψ 1 el X * + : 12 mat 0 if Δy Y max Y min - ν / (2^h 2^+ sqrt eval : Δy' Δy h 1 -^* : 12 mat x y y' fy x y Δy + y' f x y Δy - y' f - 2 Δy * / : x y y' fy' x y y' Δy' + f x y y' Δy' - f - 2 Δy' * / : A ν 1 matrix eval : D A : C A : B A : Y' A : j 2 n range : 16 mat Iter 1 dsol MaxIters Options str2num range Y' j el Y j 1 + el Y j 1 - el - 2 h * / : Y' 1 el β1 0 ≡ Y 2 el Y 1 el - h / c1 β1 / α1 β1 / Y 1 el * - if : Y' ν el β2 0 ≡ Y ν el Y n el - h / c2 β2 / α2 β2 / Y ν el * - if : 13 mat r 1 ν range P r el Q r el Φ r el 13 mat X r el Y r el Y' r el fy X r el Y r el Y' r el fy' X r el Y r el Y' r el f 13 mat eval : for A D C 13 mat 2 - h 2^P * - (2 h Q * - 2 / (2 h Q * + 2 / (13 mat eval : B j el 2 Y j el * - Y j 1 + el + Y j 1 - el + h 2^Φ j el * - eval : A 1 el D 1 el C 1 el B 1 el 14 mat β1 0 ≡ 100 Y 1 el c1 α1 / - 14 mat A 1 el 2 h * α1 β1 / * B 1 el * + (022 Y 2 el * 21 - h α1 β1 / * + (* Y 1 el * + 2 h * c1 β1 / * - h 2^Φ 1 el * - (14 mat if eval : A ν el D ν el C ν el B ν el 14 mat β2 0 ≡ 100 Y ν el c2 α2 / - 14 mat A ν el 2 h * α2 β2 / * C ν el * - (202 Y n el * 21 h α2 β2 / * + (* Y ν el * - 2 h * c2 β2 / * + h 2^Φ ν el * - (14 mat if eval : Yo A D C B - TriDiagSys : Y Y Yo + eval : 12 mat Yo norme ε Y norme * < break continue if 81 line for Iter dsol MaxIters Options str2num < X Y Y' augment eval MaxIters  rached. error if 91 line :

dsol(DE,X,tf,N) try to solve the Differential Equation system DE for the functions in X. Ev# Yv# TF# N# dsol m kg s K A cd 16 mat 11111116 mat : u# v# t# EQ_F# IG_F# EQ_F# RK_D# Clear ky# 0 : kn# 0 : kv# 0 : 14 mat ys# 0 : U# 0 : V# 0 : K# 0 : NK# 0 : IC# 03 matrix eval : 16 mat M# 0 : ES# Ev# num2str : ES'# ES#' strsplit rows : 13 mat BIC# 0 : bic# 0 : too# : 13 mat s# Yv# sys2mat.1 num2str vectorize eval a# b# 12 mat s# (findstr 1 el s#) findstr 1 el 12 mat : ts# s# a# 1 + b# a# - 1 - substr : Assuming all ≡ 12 mat ys# ky# ky# 1 + : el s# 1 a# 1 - substr : k# 1 ES'# range ES# s# s# ('(strrep : findstr 1 - ≠ M# ky# el k# : if for 41 line for s# Ev# sys2mat.1 num2str vectorize eval g# s# ≡ strsplit : g'# s# ≈ strsplit : 12 mat dsol g# 1 el 1 findrows 0 ≡ NK# NK# 1 + : s# (ts#) concat findstr 1 - ≡ a# b# 12 mat g# 1 el (findstr 1 el g# 1 el) findstr 1 el 12 mat : to# s# a# 1 + b# a# - 1 - substr : too# ≡ (too# to# ≠ (| BIC# bic# bic# 1 + : el to# : too# to# : 12 mat if g'# length 1 ≡ ico# s# 1 a# 1 - substr : S# ico#' strrep (ts#) concat strrep : f# ys# S# Find : (0 ≠ IC# IC# f# 1 el ico#' strsplit rows g# 2 el 13 mat stack : if K# NK# 1 el IC : K# NK# 2 el g# 1 el str2num g# 2 el str2num - : bvp 13 mat 41 line K# NK# 1 el IG : K# NK# 2 el g'# 1 el str2num g'# 2 el str2num - : 12 mat if 41 line g# length 1 ≡ K# NK# 1 el EQ : K# NK# 2 el s# str2num : 12 mat K# NK# 1 el EQ : K# NK# 2 el g# 1 el str2num g# 2 el str2num - : 12 mat if if 21 line dsol s# Options if 21 line for ic# 1000 IC# 1 col * IC# 2 col + IC# 3 col augment : ic# ic# 1 csort 2 col eval : ic# ic# str2num vectorize : 12 mat n# 1 ky# range kn# kn# 1 + : a# ys# n# el (ts#) concat : 12 mat K# K# a# u# kn# el Equrep : m# 1 : m# M# n# el ≤ m# m# 1 + : a# a# ('(strrep : m# M# n# el 1 - ≤ kn# kn# 1 + : kv# kv# 1 + : V# kv# el u# kn# el : 13 mat K# K# a# u# kn# el Equrep : 21 line if 21 line for kv# kv# 1 + : U# n# el v# n# el : V# kv# el v# n# el : 13 mat K# K# a# v# n# el Equrep : 51 line for EQ# K# EQ 1 findrows 2 col : v# EQ_F# EQ# eval : 12 mat bic# 1 > IC# K# IC 1 findrows 2 col : u# FI# IC# : 12 mat M# u1# u2# 21 mat FI# u1# diff vectorize u1# u2# 21 mat FI# u2# diff vectorize 0021 mat FI# - augment eval : dsol dsolver Options lbvp.2 ≡ u# v# FE# EQ# 1 el : dv# u1# u2# 21 mat v1# 11 mat FE# v1# diff : 12 mat L# u1# u2# 21 mat v1# 11 mat FE# u2# diff dv# / u1# u2# 21 mat v1# 11 mat FE# u1# diff dv# / 0021 mat 011 mat FE# dv# / - 13 mat : lf#(ts#): L# num2str concat str2num v1# lf# BIC# str2num vectorize M# N# lbvp.2 41 line dsol dsolver Options bvp.2 ≡ vg# dsol DefaultGuess Options str2num : vg# vg# vg# stack : f# K# IG 1 findrows : num2str 0 ≠ IG# f# 2 col : u# IG_F# IG# eval : 12 mat vg# vg# 01020 -^IG_F# al_nleqsolve 2 el N# 1 + 1 matrix + eval : 21 line vg# 0 : if v# BF# EQ# : ε# dsol eps Options str2num : go# dsol DefaultGuess Options str2num 11 mat : t# y# y'# bf# ts# : t# num2str concat str2num u# y# y'# 12 mat : go# BF# al_nleqsolve 1 el 31 line : bf# BIC# str2num vectorize vg# M# N# ε# bvp.2 81 line Not Implemented. error if 31 line vg# dsol DefaultGuess Options str2num ky# 1 matrix + eval : f# K# IG 1 findrows : num2str 0 ≠ IG# f# 2 col : v# IG_F# IG# eval : vg# vg# 01020 -^IG_F# al_nleqsolve : 13 mat if v# vg# : dsol UpdateGuess Options 1 ≡ t# u# v# RK_D# ts# : t# num2str concat str2num v# v# EQ_F# al_nleqsolve : V# eval 31 line : t# u# RK_D# ts# : t# num2str concat str2num v# vg# EQ_F# al_nleqsolve : V# eval 31 line : if dsol dsolver Options (ic#,str2num(to#),TF#,N#,RK_D#) concat str2num 51 line if 121 line :

DE# Y# N# dsol DE# Y# 0 N# dsol 11 line :

dsol dsolver rkfixed steps 100 nsolver al_nleqsolve stepsmax 0 eps 10^{-7} DefaultGuess 1 UpdateGuess 1 MaxIters 5082 mat :

nsol nsolver al_nleqsolve stepsmax 0 eps 10^{-7} 32 mat :

dsol(DE,X,tf,N) try to solve the Differential Equation system DE for the functions in X.

N 100 :

<content> Example </content>

Example

Default solver is

dsol dsolver Options rkfixed

It could be changed with

dsol dsolver solver ≡ Options

or in the solver block

t us 0.5 t * t 2^1 + sqrt + : t vs 2 t * 4 t * cos * : 21 sys

RK t u t v 2 1 sys 4 N dsol : RK 1 col RK 2 col augment x us 2 1 sys

0 u 1 ≡

0 v 0 ≡

0 v' 2 ≡

dsolver Adams ≡

t v'' t u t v' abs sqrt * / - 322 t * 4 t * cos * 4 t * sin + (* t 21 t 2^+ sqrt * + (24 t * 4 t * sin * 4 t * cos - (* abs sqrt * / ≡

t u' t v' * t u / - 21 t 2^+ sqrt 2 t * + (* 4 t * 4 t * sin * 4 t * cos - (* 1 t 2^+ sqrt t 21 t 2^+ sqrt * + (* / ≡

RK 1 col RK 3 col augment x vs 2 1 sys <content> Example (Exact ODE) </content>

Example

ε 0.01 :

x ys 1 x x ln * / :

yo 1 ε + ys :

RK x y 6 N dsol : RK 1 col RK 2 col augment x ys 2 1 sys

x y x ln * x y + x x ln * x y' * + 0 ≡

1 ε + y yo ≡

Example

C 2 - :

x ys x 2^5 x * y * + y 3^+ C ≡ y 1 - roots 11 line :

yo 0 ys :

RK x y 8 N dsol : RK 1 col RK 2 col augment Xs RK 1 1 el RK RK rows 1 el 10 pR : Xs ys vectorize . augment 2 1 sys

2 x * 5 x y * + 3 x y (2^* 5 x * + (x y' * + 0 ≡

0 y yo ≡

<content> Example Guess for the higest derivatives </content>

Example

This example shows how to use the guess value for the higest derivatives

x ys 19 / x 3^9 + 63 x 3^* sqrt * + (* :

RK x y 3 N dsol : RK 1 col RK 2 col augment x ys 2 1 sys

x y' 2^x x y * - 2 x * ≡

dsolver Rkadapt ≡

0 y 1 ≡

Above solve block use

dsol DefaultGuess Options 1

It could be changed with

x ys 19 / x 3^9 + 63 x 3^* sqrt * - (* :

dsol DefaultGuess 1 - ≡ Options

or specifiying a value in the solve block

RK x y 3 N dsol : RK 1 col RK 2 col augment x ys 2 1 sys

x y' 2^x x y * - 2 x * ≡

dsolver Rkadapt ≡

0 y 1 ≡

0 y' 1 - ≈

<content> Example </content>

Example

xo 1 :

yo 2 :

x C eq 2 x * e x 2^^* e x 2^^C + / ln :

Co xo C eq yo - C 2 - roots :

x ys x Co eq :

RK x y 10 N dsol : RK 1 col RK 2 col augment x ys 2 1 sys

x y' e x y^+ 2 x * 1 x / + ≡

xo y yo ≡

<content> Example </content>

Example

α 12 :

σ us 1 α σ * + ln : σ vs e σ -^1 + : 21 sys

RK σ u σ v 2 1 sys 4 N dsol : RK 1 col RK 2 col augment x us 2 1 sys

σ v'' σ v * σ σ u * - e σ^1 + σ 1 α σ * + ln * e 2 σ *^* - e 2 σ *^/ ≡

σ u' σ v'' * 2 σ * cos + α 2 σ * cos 1 α σ * + (* e σ^* + 1 α σ * + (e σ^* / ≡

0 v 2 ≡

0 v' 1 - ≡

0 v'' 1 ≈

RK 1 col RK 3 col augment x vs 2 1 sys

0 u 0 ≡

0 u' α ≈

dsolver dn_AdamsMoulton ≡

<content> Example Abel </content>

Example

x ys 1 x y * - (e y y 2 x 3^* + 2 x 2^* / *^* 0 ≡ y 1 roots 11 line :

RK x y 8 N dsol : RK 1 col RK 2 col augment Xs RK 1 1 el RK RK rows 1 el 10 pR : Xs ys vectorize . augment 2 1 sys

x x x y * x 4^+ 1 - (* x y' * x y x x y * x 4^- 1 - (* ≡

1 y 1 ≡

1 y' 1 ≈

dsolver Rkadapt ≡

<content> Example </content>

Example

t xs 2 t sin * : t ys t e t -^* : 21 sys

RK t x t y 2 1 sys 10 N dsol : RK 1 col RK 2 col augment x xs 2 1 sys

0 x 0 ≡

0 x' 2 ≡

0 y 0 ≡

0 y' 1 ≡

t y'' t * t x t y * - t 2 t - 2 t sin * + (* e t^/ +

t x'' 4 t * cos t y' * - 2 t sin * e t^* 4 t * cos 1 t - (* + e t^/ +

RK 1 col RK 4 col augment x ys 2 1 sys

0 x'' 0 ≈

0 y'' 2 - ≈

<content> Example </content>

Example

t xs t 2^6 - : t ys e t sin^: t zs 0.5 t 2^* 2 - : 31 sys

RK t a t b t c 3 1 sys 10 N dsol : RK 1 col RK 2 col augment x xs 2 1 sys

0 a 6 - ≡

0 b 1 ≡

0 c 2 - ≡

0 a' 0 ≡

0 b' 1 ≡

0 c' 0 ≡

t a'' t b * t c - 41 e t sin^+ (* t 2^- 2 / ≡

t b'' t t c' * - e t sin^t cos (2^t sin - (* t 2^- ≡

RK 1 col RK 4 col augment x ys 2 1 sys

t c'' t a'' * 3 t b * - 23 e t sin^* - ≡

0 a'' 2 ≈

0 b'' 1 ≈

0 c'' 1 ≈

dsolver Rkadapt ≡

RK 1 col RK 6 col augment x zs 2 1 sys <content> Example </content>

Example

x ys 11 Yv 01 Yv - (1 x Jv * 11 Jv 01 Jv - (1 x Yv * - 11 Yv 01 Jv * 01 Yv 11 Jv * - / :

RK x y 7 5 N * dsol : RK 1 col RK 2 col augment x ys 2 1 sys

x 2^x y''' * x 3^- 3 x * + (x y'' * + x 3^- 3 x * + (x y * + 0 ≡

1 y 1 ≡

1 y' 0 ≡

1 y'' 0 ≡

<content> Example UpdateGuess </content>

The option UpdateGuess = true is used for update the value of the y' guess at each iteration. Default value is true (1). First example requires it set to false, and second to true.

Example

dsol UpdateGuess 0 ≡ Options 0

t ys 0.5 π * 1 e t -^- (* :

RK t y 4 N dsol : RK 1 col RK 2 col augment x ys 2 1 sys

t y' sin t y cos ≡

0 y 0 ≡

0 y' 0.5 ≈

Example

dsol UpdateGuess 1 ≡ Options 1

C 0 :

t ys y 2^1 + sqrt 12 / y 2^1 + sqrt 1 - y 2^1 + sqrt 1 + / ln * + t C + ≡ y 1 roots 11 line :

yo 0 ys : 0.6627

RK t y 4 N dsol : RK 1 col RK 2 col augment Xs RK 1 1 el RK RK rows 1 el 10 pR : Xs ys vectorize . augment 2 1 sys

t y' asin t y atan ≡

0 y yo ≡

0 y' 0 ≈

Boundary Values Problem

lbvp.2

x f ab M N lbvp.2

solves the linear ODE

ab a b 12 mat ≡

y'' p y' * + q y * + r ≡

subject to the
boundary conditions

α 1 el a y * β 1 el a y' * + c 1 el ≡ α 2 el b y * β 2 el b y' * + c 2 el ≡ 21 sys

f is such that

x f x p x q x r 13 mat ≡

and

M α1 β1 c1 α2 β2 c2 23 mat ≡

bvp.2

x y y' φ x Yo M N ε bvp.2

solves the non linear ODE

x a b 12 mat ≡

y'' x y y' φ ≡

subject to the the same boundary conditions, with Yo as guess for the solution with dimension N+1. If Yo≡0, bvp.2 try with a line between f(a) and f(b) as guess.
ε is used as the tolerance for a Newton solver.

<content> Short hands for plots </content>

Short hands for plots

x xb X Y Y' x y x y' G Ys Y's Clear 1

Err X Y Y' 13 mat sol Cols : Y ψ X feval vectorize - norme 21 line :

x ψ' x ψ x diff 11 line :

str sol Plot X Y Y' 13 mat sol 1 col sol 2 col sol 3 col 13 mat : Xs X 1 el X X rows el 300 pR : str XY ≡ u v us vs 14 mat X Y Xs ψ Xs feval vectorize 14 mat eval : u v us vs 14 mat Y Y' ψ Xs feval vectorize ψ' Xs feval vectorize 14 mat eval : if us vs augment u 1 el u u rows el 12 mat v 1 el v v rows el 12 mat pStem u v . 8 augment 31 sys 41 line :

<content> lbvp example </content>

Example

a b 12 mat 0512 mat :

h x x cos * :

h' h x diff :

Numeric
solution

sol x y N dsol :

x y'' 2 h * x y' * + h 2^h' + (x y * + 0 ≡

a y 3 a y' * + 0.1 ≡

b y' 6 - ≡

dsolver lbvp.2 ≡

Analytic
solution

x ψ A B x * + (e x x sin * x cos + (-^* :

A B 21 mat a ψ 3 a ψ' * + 0.1 ≡ b ψ' 6 - ≡ 21 mat A B 21 mat roots :

Err 0.16

XY sol Plot YY' sol Plot <content> lbvp example </content>

Example

sol x y N dsol :

x y'' 3 x y' * + 4 x y * - 0 ≡

dsolver lbvp.2 ≡

Numeric
solution

0 y 1 ≡

1 y 1 ≡

Analytic
solution

x ψ e 44 x * - (^e x^+ e x 1 + (^+ e x 2 + (^+ e x 3 + (^+ (1 e + e 2^+ e 3^+ e 4^+ (/ :

Err 5.43 105 -^*

XY sol Plot YY' sol Plot <content> lbvp example </content>

Example

sol t x N dsol :

t x'' 3 t x' * + 6 t x * + 5 ≡

1 x 3 ≡

Numeric
solution

dsolver lbvp.2 ≡

3 x 23 x' * + 5 ≡

Analytic
solution

t ψ 0.0696737 e 1.5 t * -^* 11.9605 e 1.5 t *^* 1243.5 1.93649 t * sin * + 2857.67 1.93649 t * cos * + (* :

Err 0.21

XY sol Plot YY' sol Plot <content> lbvp example </content>

Example

p 105 -^:

a b 12 mat 0.1 - 0.1 12 mat :

Analytic
solution

x ψ x p x 2^+ sqrt / :

Numeric
solution

sol x y 1000 dsol :

Err 0.1

x y'' 3 p * p x 2^+ (2^/ x y * + 0 ≡

a y a ψ ≡

dsolver lbvp.2 ≡

b y b ψ ≡

XY sol Plot YY' sol Plot <content> lbvp example </content>

Example

a b 12 mat 13 π * / 3 π / 12 mat eval :

Numeric
solution

sol x y N dsol :

dsolver lbvp.2 ≡

x y'' 2 x / x y' * + x y x 4^/ + 0 ≡

a y 0 ≡

b y 3 sqrt 2 / eval ≡

Analytic
solution

x ψ 1 x / sin :

Err 0.32

XY sol Plot YY' sol Plot <content> bvp with two solutions </content>

Example

a b 12 mat 0412 mat :

y.b 2 - :

Numeric
solution

sol x y N dsol :

X Y Y' 13 mat sol Cols :

x y'' x y abs + 0 ≡

a y 0 ≡

dsolver bvp.2 ≡

y'.a Clear 1

y'.a 0 :

Seems like a bug in the component.

b y y.b ≡

This problem have two solutions. I compare the bvp solution with the shooting method, because can't found a symbolic expression.

Shoot
solution

y'.a shoot x ψ b N Rkadapt :

x ψ'' x ψ abs + 0 ≡

a ψ 0 ≡

a ψ' y'.a ≡

y'.a Eq y'.a shoot N 1 + 2 el y.b - :

y'.a Y' 1 el Eq al_nleqsolve 1 el : 0.0733 -

Ξ Ψ Ψ' 13 mat y'.a shoot Cols :

x ψ' Ξ Ψ' x cinterp :

Err 0

XY sol Plot YY' sol Plot sol x y N dsol :

For the other solution: call bvp with a new guess

x y'' x y abs + 0 ≡

a y 0 ≡

dsolver bvp.2 ≡

X Y Y' 13 mat sol Cols :

0 y' 1 ≈

b y y.b ≡

Solve the shoot with a new guess from bvp

y'.a Y' 1 el Eq al_nleqsolve 1 el : 2.0666

Ξ Ψ Ψ' 13 mat y'.a shoot Cols :

x ψ' Ξ Ψ' x cinterp :

Err 0

XY sol Plot YY' sol Plot

Alvaro

4 appVersion 1.0.8348.30405