Re: A problem with ilaplace()

[Leo Baumann]

Am 14.11.2017 um 18:22 schrieb Doug Stewart:

On Tue, Nov 14, 2017 at 9:47 AM, Leo Baumann <address@hidden> wrote:

Hello folks,

I have reduced the problem. -

It is explained in Octave script:

#--------------------------------------------------------------------------

clc;
clear all;
pkg load symbolic;
syms p t;

#a few simple calculations
Rgen=50;C1=100e-9;Rlast=50;
fg=10e3;
R1=1/(2*pi*10e3*C1)-Rgen;
wg=2*pi*fg;

#this works fine, Tp1 is a simple imagined example ...
Tp1=1/(5*p+6);
laplace1=1/p*Tp1;
ux=function_handle(ilaplace(laplace1, p, t));

#this doesn't work, I don't understand that because Tp2 is the same
#as the example with other numbers
Tp2=1/(1+(Rgen+R1)/Rlast+p*C1*(Rgen+R1));
laplace2=1/p*Tp2;
ux=function_handle(ilaplace(laplace2, p, t));

#----------------------------------------------------------------------------------------

I don't understand that the second example doesn't work, all variables are defined.

Any hint?

Thanks - Leo


_

What results are you looking for. What is your overall goal?

At the moment my script makes no use. I was preparing a general usable script to evaluate step function responses for electric networks. I could do that numeric with following Octave function, but the symbolic answer is much more interesting. -

# Talbot suggested that the Bromwich line be deformed into a contour that begins
# and ends in the left half plane, i.e., z ? -8 at both ends.
# Due to the exponential factor the integrand decays rapidly
# on such a contour. In such situations the trapezoidal rule converge
# extraordinarily rapidly.
# For example here we compute the inverse transform of F(s) = 1/(s+1) at t = 1
#-----------------------------------------------------------------------------
#>> pkg load symbolic
#>> syms s
#>> F=1/(s+1)
#F = (sym)
#
#    1
#  -----
#  s + 1
#
#>> error=talbot(function_handle(F),1,24)-exp(-1)
#ans =   1.6098e-015
#-------------------------------------------------------------------------------
#
# Talbot method is very powerful here we see an error of 1.61e-015
# with only 24 function evaluations
#
# Created by Fernando Damian Nieuwveldt     
# email:address@hidden
# Date : 25 October 2009
#
# Reference
# L.N.Trefethen, J.A.C.Weideman, and T.Schmelzer. Talbot quadratures
# and rational approximations. BIT. Numerical Mathematics,
# 46(3):653 670, 2006.

function ilt = talbot(f_s, t, N)
h=2*pi/N;

#   Shift contour to the right in case there is a pole on the positive real axis : Note the contour will
#   not be optimal since it was originally devoloped for function with
#   singularities on the negative real axis
#   For example take F(s) = 1/(s-1), it has a pole at s = 1, the contour needs to be shifted with one
#   unit, i.e shift  = 1. But in the test example no shifting is necessary
 shift=0;

 ans=0;
 for k=0:N
   theta=-pi+(k+1/2)*h;
   z=shift+N/t*(0.5017*theta*cot(0.6407*theta)-0.6122+0.2645i*theta);
   dz=N/t*(-0.5017*0.6407*theta/sin(0.6407*theta)^2+0.5017*cot(0.6407*theta)+0.2645i);
   ans=ans+exp(z*t)*f_s(z)*dz;
 endfor
   ilt=real((h/(2i*pi))*ans);
endfunction

Regards - Leo