How to define Recurrence Relations in IGOR
shankar.dutt
Hi,
I don't have much experience with programming, but I have recurrence relations for which there is no closed-form. I want to define them in IGOR Pro. I used simple if..else to define them but, when I try to get an output I get the following error:
"out of stack space (too much macro or function recursion)". It just works for P0(0), P0(1), Q0(0) and Q0(1). What would be the way to go?
Please help.
Function P0(n)
variable n
variable a,q
if(n==0)
return 0
elseif(n==1)
return 1
elseif(n>=2)
return (2*P0(n-1)-2*a*P0(n)-q*Q0(n))/(n+1)
endif
End
Function Q0(n)
variable n
variable a,q
if(n==0)
return 0
elseif(n==1)
return 0
elseif(n>=2)
return (2*Q0(n-1)-2*a*Q0(n)-q*P0(n))/(n)
endif
End
variable n
variable a,q
if(n==0)
return 0
elseif(n==1)
return 1
elseif(n>=2)
return (2*P0(n-1)-2*a*P0(n)-q*Q0(n))/(n+1)
endif
End
Function Q0(n)
variable n
variable a,q
if(n==0)
return 0
elseif(n==1)
return 0
elseif(n>=2)
return (2*Q0(n-1)-2*a*Q0(n)-q*P0(n))/(n)
endif
End
In your example above I can't see where the variables a and q are assigned any value?
Note also, q has a special meaning in Igor - it is the column index of a destination wave when used in a multidimensional wave assignment. See DisplayHelpTopic("Waveform Arithmetic and Assignment").
Both of your functions are (in the case of n >= 2) functions of themselves. i.e. the function calls itself *with the same parameter value*. This leads to an infinite nest of function calls. For example, calling P0(2) will (in the last return statement) calls a new instance of the function P0(2), which in turn calls a new instance of P0(2), and so on.
I have no idea whether the following is valid, but...
if you take the function assignment P0(n) = (2*P0(n-1)-2*a*P0(n)-q*Q0(n))/(n+1),
you can rearrange this to something like:
P0(n) = (2*P0(n-1) - q*Q0(n)) / (n + 1 + 2*a)
Does this help?
Regards,
Kurt
September 19, 2019 at 01:47 am - Permalink
OP approach is less than ideal in that there are limits to recurrence (as OP found) and it is very inefficient because one is unnecessarily calculating previous values multiple times. One approach is to estimate the largest parameter N that one wants to support and then create a couple of waves say:
Make/n=(N) pWave=nan,qWave=nan
From here define your functions similar to:
variable n
Wave pWave
variable a,q,res
if(numType(pWave[n])!=2)
return pWave[n]
endif
if(n==0)
pWave[0]=0
return 0
elseif(n==1)
pWave[1]=1
return 1
elseif(n>=2)
res=(2*P0(n-1,pWave)-2*a*P0(n,pWave)-q*Q0(n,pWave))/(n+1)
pWave[n]=res
return res
endif
End
I hope this helps,
A.G.
September 19, 2019 at 06:47 am - Permalink