Archive for November, 2017

Convincing maple to solve an ODE with Neumann conditions at a symbolic valued location

Friday, November 17th, 2017

I can use maple to solve a 1D second-order ODE with Dirichlet boundary conditions at symbolic-valued locations:

# Z'' = 0, Z(a)=0, Z(b) = 1
dsolve({diff(Z(r),r,r) = 0,Z(a)=0,Z(b)=1}); 

This correctly returns

                  r       a
       Z(r) = - ----- + -----
                a - b   a - b

I can also easily convince maple to solve this ODE with some Neumann (normal derivative) boundary conditions at at fixed-value, numeric location:

# Z'' = 0, Z(a) = 1, Z'(0) = 0
dsolve({diff(Z(r),r,r) = 0,Z(a)=1,eval(diff(Z(r),r),r=0)=0});


                 Z(r) = 1

But if I try naively to use a Neumann condition at a symbolic value location

# Z'' = 0, Z(a) = 1, Z'(b) = 0
dsolve({diff(Z(r),r,r) = 0,Z(a)=1,eval(diff(Z(r),r),r=b)=0});

then I get an error:

Error, (in dsolve) found differentiated functions with same name but depending on different arguments in the given DE system: {Z(b), Z(r)}

After a long hunt, I found the solution. dsolve takes an optional second argument that can tell it what the dependent variable actually is. So the correct call is:

# Z'' = 0, Z(a) = 1, Z'(b) = 0
dsolve({diff(Z(r),r,r) = 0,Z(a)=1,eval(diff(Z(r),r),r=b)=0});

and this gives the correct answer

                 Z(r) = 1

MATLAB gotcha inverting a (sparse) diagonal matrix

Thursday, November 2nd, 2017

Just got burned by a silly Matlab gotcha. Suppose you have a diagonal matrix D and you want to compute the inverse square root matrix:

Disqrt = diag(1./sqrt(diag(D))

But this will be dense!

Instead use

Disqrt = diag(sqrt(diag(D).^-1)

Or maybe

Disqrt = diag(diag(D).^-0.5)

Not sure if there’s an accuracy difference (hopefully not).