help-octave
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: complex line integral in octave


From: Stephen Montgomery-Smith
Subject: Re: complex line integral in octave
Date: Mon, 15 Apr 2013 17:08:48 -0500
User-agent: Mozilla/5.0 (X11; Linux i686; rv:17.0) Gecko/20130329 Thunderbird/17.0.5

On 04/15/2013 05:02 PM, Stephen Montgomery-Smith wrote:
> This is a programming issue, not mathematics.  If you look at "help
> quad" you will see that the function takes only one parameter.  The
> other parameters, La,Lb,Lc,Ia,Ib,Ic,RC, will have to be passed in some
> other manner.  My suggestion would be to use global variables.  But
> others more proficient in octave or matlab may have better suggestions.
> 
> Also, if you want to use other integration algorithms, some of them
> require the function to accept vectors, and then evaluate the function
> at all those values.  And you may as well try them out to see if they
> give better or faster results.
> 
> Also octave can handle complex numbers, so I wouldn't mess with
> splitting up the real and imaginary parts.
> 
> To compute the line integral from z=a to z=b, I would use the following
> code.  But it isn't tested.  I also pass a and b using global variables.
> 
> function g = myf(t)
> global La,Lb,Lc,Ia,Ib,Ic,RC
> global a,b
> 
> a0=koeff0(La,Lb,Lc,Ia,Ib,Ic,RC);
> a1=koeff1(La,Lb,Lc,Ia,Ib,Ic,RC);
> a2=koeff2(La,Lb,Lc,Ia,Ib,Ic,RC);
> a3=koeff3(La,Lb,Lc,Ia,Ib,Ic,RC);
> a4=koeff4(La,Lb,Lc,Ia,Ib,Ic,RC);
> a5=koeff5(La,Lb,Lc,Ia,Ib,Ic,RC);
> a6=koeff6(La,Lb,Lc,Ia,Ib,Ic,RC);
> 
> z = a + (b-a)*t;
> 
> f = (((((a6.*z.+a5).*z.+a4).+a3).*z.+a2).*z.+a1).*z.+a0;
> fp = ((((6*a6.*z.+5*a5).*z.+4*a4).+3*a3).*z.+2*a2).*z.+a1;
> 
> g = z.*fp./f

I mean

g = z.*fp./f.*(b-a)

> 
> end
> 
> On 04/15/2013 08:51 AM, Urs Hackstein wrote:
>> Thanks a lot for your explanations, they were very helpful. I have still
>> a problem concerning the numerical integration. I want to do the
>> computation for several different polynomials f, so I have seven
>> parameters La,Lb,Lc,Ia,Ib,Ic, RC besides the integration variable t. My
>> attempt for the integral over the first edge of the rectangle was the
>> following:
>>
>> function g1=complexfunction1(La,Lb,Lc,Ia,Ib,Ic,RC,t)
>> sd=z(1,-10.^11,1,poly(La,Lb,Lc,Ia,Ib,Ic,RC));
>> wmax=3*10.^11;
>> a0=koeff0(La,Lb,Lc,Ia,Ib,Ic,RC);
>> a1=koeff1(La,Lb,Lc,Ia,Ib,Ic,RC);
>> a2=koeff2(La,Lb,Lc,Ia,Ib,Ic,RC);
>> a3=koeff3(La,Lb,Lc,Ia,Ib,Ic,RC);
>> a4=koeff4(La,Lb,Lc,Ia,Ib,Ic,RC);
>> a5=koeff5(La,Lb,Lc,Ia,Ib,Ic,RC);
>> a6=koeff6(La,Lb,Lc,Ia,Ib,Ic,RC);
>> G1=ones(1,6).*(sd-10.^8+2*10.^8*t).^(1:6);
>> H1=[a1 2*a2 3*a3 4*a4 5*a5 6*a6];
>> G11=ones(1,7).*(sd-10.^8+2*10.^8*t).^(0:6);
>> H11=[a0 a1 a2 a3 a4 a5 a6]
>> g1=(7*(sd-10.^8+2*10.^8*t).^7+H1.*G1')*2*10.^8/((sd-10.^8+2*10.^8*t).^7+H11*G11')
>> endfunction
>>
>>  Here sd is the real part of the root inside the rectangle and the ai
>> are the coefficients of the polynomial.
>>
>> function g1r=real_part1(La,Lb,Lc,Ia,Ib,Ic,RC,t)
>> g1r=real(complexfunction1(La,Lb,Lc,Ia,Ib,IC,RC,t))
>> endfunction
>>
>> function g1i=imag_part1(La,Lb,Lc,Ia,Ib,Ic,RC,t)
>> g1i=imag(complexfunction1(La,Lb,Lc,Ia,Ib,IC,RC,t))
>> endfunction
>>
>> function complex_answer1=complex1(La,Lb,Lc,Ia,Ib,Ic,RC)
>> real_integral1=quad("real_part1",0,1)
>> imag_integral1=quad("imag_part1",0,1)
>> complex_answer1=real_integral1+I*imag_integral1
>> endfunction
>>
>> I guess that quad doesn't know which variable is the integration
>> variable, but how can we fix it?
>> Thanks a lot in advance.
>>
>>
>> 2013/4/12 Stephen Montgomery-Smith <address@hidden
>> <mailto:address@hidden>>
>>
>>     On 04/12/2013 05:29 AM, Urs Hackstein wrote:
>>     > Dear Stephen,
>>     >
>>     > thanks a lot for your long reply. Unfortunately, we are talking about
>>     > different integrals. Your are dealing with the zero-counting integral
>>     > \int_\alpha f ' (s)/f(s) ds, I work on the integral \int_\alpha s*f '
>>     > (s)/f(s) ds, where * is the common multiplication. From complex
>>     function
>>     > theory we know that it equals a, if a is the only root of f in the
>>     > rectangle. We want to compute this root, thus method 2 doesn't
>>     work for
>>     > us. Remains method 1: what is the best (most efficient and stable)
>>     built
>>     > in integration command?
>>     >
>>
>>     First, you meant 2 pi I a, not a.  I'm sure it was a typo, but it is the
>>     kind of typo that propagates into code and creates a hard to trace bug.
>>      (At least, that is what it does for me.)
>>
>>     Second, I think you have a fairly straightforward integral to compute,
>>     so I think any numerical method should do fine.  (Unless one of the
>>     roots lies very close or is on the curve.)
>>
>>     Third, if you know that there is only one root inside the rectangle, you
>>     could use the output of Method 1 as the input to Newton's Method, which
>>     will very quickly and accurately converge to the desired root.
>>
>>     Finally, friends of mine highly recommend the NIntegrate command in
>>     recent versions of Mathematica (version 8 or higher).  It might even
>>     handle the case where one root lies on the curve (where presumably it
>>     correctly computes the principle value).  So perhaps you could check
>>     your answers using the wolfram alpha web site,
>>
>>
>>
>> 2013/4/12 Stephen Montgomery-Smith <address@hidden
>> <mailto:address@hidden>>
>>
>>     On 04/12/2013 05:29 AM, Urs Hackstein wrote:
>>     > Dear Stephen,
>>     >
>>     > thanks a lot for your long reply. Unfortunately, we are talking about
>>     > different integrals. Your are dealing with the zero-counting integral
>>     > \int_\alpha f ' (s)/f(s) ds, I work on the integral \int_\alpha s*f '
>>     > (s)/f(s) ds, where * is the common multiplication. From complex
>>     function
>>     > theory we know that it equals a, if a is the only root of f in the
>>     > rectangle. We want to compute this root, thus method 2 doesn't
>>     work for
>>     > us. Remains method 1: what is the best (most efficient and stable)
>>     built
>>     > in integration command?
>>     >
>>
>>     First, you meant 2 pi I a, not a.  I'm sure it was a typo, but it is the
>>     kind of typo that propagates into code and creates a hard to trace bug.
>>      (At least, that is what it does for me.)
>>
>>     Second, I think you have a fairly straightforward integral to compute,
>>     so I think any numerical method should do fine.  (Unless one of the
>>     roots lies very close or is on the curve.)
>>
>>     Third, if you know that there is only one root inside the rectangle, you
>>     could use the output of Method 1 as the input to Newton's Method, which
>>     will very quickly and accurately converge to the desired root.
>>
>>     Finally, friends of mine highly recommend the NIntegrate command in
>>     recent versions of Mathematica (version 8 or higher).  It might even
>>     handle the case where one root lies on the curve (where presumably it
>>     correctly computes the principle value).  So perhaps you could check
>>     your answers using the wolfram alpha web site,
>>
>>
>>
>>
>> 2013/4/12 Stephen Montgomery-Smith <address@hidden
>> <mailto:address@hidden>>
>>
>>     On 04/12/2013 05:29 AM, Urs Hackstein wrote:
>>     > Dear Stephen,
>>     >
>>     > thanks a lot for your long reply. Unfortunately, we are talking about
>>     > different integrals. Your are dealing with the zero-counting integral
>>     > \int_\alpha f ' (s)/f(s) ds, I work on the integral \int_\alpha s*f '
>>     > (s)/f(s) ds, where * is the common multiplication. From complex
>>     function
>>     > theory we know that it equals a, if a is the only root of f in the
>>     > rectangle. We want to compute this root, thus method 2 doesn't
>>     work for
>>     > us. Remains method 1: what is the best (most efficient and stable)
>>     built
>>     > in integration command?
>>     >
>>
>>     First, you meant 2 pi I a, not a.  I'm sure it was a typo, but it is the
>>     kind of typo that propagates into code and creates a hard to trace bug.
>>      (At least, that is what it does for me.)
>>
>>     Second, I think you have a fairly straightforward integral to compute,
>>     so I think any numerical method should do fine.  (Unless one of the
>>     roots lies very close or is on the curve.)
>>
>>     Third, if you know that there is only one root inside the rectangle, you
>>     could use the output of Method 1 as the input to Newton's Method, which
>>     will very quickly and accurately converge to the desired root.
>>
>>     Finally, friends of mine highly recommend the NIntegrate command in
>>     recent versions of Mathematica (version 8 or higher).  It might even
>>     handle the case where one root lies on the curve (where presumably it
>>     correctly computes the principle value).  So perhaps you could check
>>     your answers using the wolfram alpha web site,
>>
>>
> 
> _______________________________________________
> Help-octave mailing list
> address@hidden
> https://mailman.cae.wisc.edu/listinfo/help-octave
> 
> 



reply via email to

[Prev in Thread] Current Thread [Next in Thread]