Solving an algebraic-integral equation in Mathematica

  • #1

Main Question or Discussion Point

Hi,

I'm trying to solve an algebraic-integral equation (I don't know if this is the proper description of this class of equations, it just seems like the least wrong way to describe them) and have run into several issues that I'll describe below, but first I'll outline the problem.

I'm trying to solve for [itex]k[/itex] dependent variables where the governing equation for the [itex]k-th[/itex] mode takes the form

[tex] \alpha_{k} +\frac{1}{\lambda}\gamma \int_0^{\lambda} \frac{(1+\mathcal{X})\cos k\xi +\mathcal{Y}\sin k\xi }{\sqrt{(1+\mathcal{X})^2+\mathcal{Y}^2} } d\xi =0 [/tex]

where

[tex] \mathcal{X}=2\sum_{n=1}^N n(\alpha_n \cos n\xi);\quad \& \quad \quad \mathcal{Y} = 2\sum_{n=1}^N n(-\alpha_n \sin n\xi ).
[/tex]

and the parameters [itex]\lambda, \gamma [/itex] and N are prescribed, where N is a parameter that will be be increased to meet convergence criterion. The independent variable of the system is [itex]\gamma[/itex] and we would like to solve for the [itex]\alpha_k[/itex].

Now, I don't really know how to deal with the integral term and have no reason to suspect there's a better way to rewrite it, so I am doing the most naive thing that comes to mind, namely, taking the so called trapezoidal approximation of the integral over even spacing. That is, given a function [itex]f[/itex] on an interval [itex](a,b)[/itex], we have

[tex]\int_a^b f(x) dx =\frac{b-a}{2M} \left(f(x_1)+2\sum_{i=2}^Mf(x_i) +f(x_{M+1})\right) +R_M[/tex]

where [itex]x_i=a+\frac{b-a}{M}*i [/itex].

My goal is to increase the resolution of the integral until the remainder [itex] R_M[/itex] is sufficiently small.

I'm having a few problems. My main question is, is this the proper way to go about problems of this nature? Something is wrong with either my physics (which I understand you cannot help with without more details, which I can provide if interested) or my code as I get very strange non-convergent behavior as M is increase. Secondly, I am finding that for M<N, the gross majority of the energy tends to congregate in the M-th mode. I'm not completely sure why the system would bias this mode so I've just been considering cases where M>N. Last, but certainly not least, this code ends up being very memory intensive for test examples that I'm playing with, so I'm looking at ways of getting around this.

Any suggestions for my code, or even on how to approach these types of problems, would be much appreciated.


My mathematica code is included below.
M = 20;
\[Kappa] = 10^6;
\[Lambda] = 2 \[Pi]/\[Kappa];
\[Gamma] = 1.0126;
f1a[k_, x_] := ((1 +
2 Sum[n*(\[Alpha][n]*Cos[2 \[Pi]*n*x]), {n, 1, M}]) Cos[
2 \[Pi]*k*x] +
2*Sin[2 \[Pi]*k*x]*
Sum[n*(-\[Alpha][n]*Sin[2 \[Pi]*n*x]), {n, 1, M}])/
Sqrt[(1 + 2 Sum[n*(\[Alpha][n]*Cos[2 \[Pi]*n*x]), {n, 1, M}])^2 +
Sum[2 n*(-\[Alpha][n]*Sin[2 \[Pi]*n*x]), {n, 1, M}]^2];
MMax = 22;
I1a[j_] := 1/(MMax) (Sum[f1a[j, (i - 1)/(MMax)], {i, 1, MMax}]);
Governing[j_] := \[Alpha][j] - \[Gamma]*I1a[j];
IC = Join[Table[{\[Alpha][1], .4}, {i, 1}],
Table[{\[Alpha], 0}, {i, 2, M}]];
S =
FindRoot[Table[With[{i = i}, Governing == 0], {i, 1, M}],IC];
Y[x_] := Evaluate[
Sum[-2*(\[Alpha][k]^2)/\[Kappa] +
2*Cos[2 \[Pi]*x*k]*\[Alpha][k]/\[Kappa], {k, 1, M}] /.
S];
X2[x_] :=
x*\[Lambda] +
Evaluate[
Sum[Sin[2 \[Pi]*x*k]*\[Alpha][k]/\[Kappa], {k, 1, M}] /. S];
ParametricPlot[{X2[x], Y[x]}, {x, 0, 1}, AspectRatio -> 1/4]


Nick
 

Answers and Replies

  • #2
SteamKing
Staff Emeritus
Science Advisor
Homework Helper
12,796
1,667
Integral equations can be tricky to solve, and the numerical methods which can solve them involve using a lot of complicated mathematics to figure out the solutions.

Because of the nature of the integral expression, the use of the trapezoidal method to calculate the integral may not provide accurate results. There are numerical integration methods which are specially adapted to give accurate results when trigonometric or other periodic functions are involved.
 

Related Threads on Solving an algebraic-integral equation in Mathematica

Replies
4
Views
3K
Replies
1
Views
716
Replies
20
Views
14K
Replies
3
Views
2K
Replies
5
Views
2K
Replies
1
Views
649
Replies
4
Views
3K
Replies
6
Views
10K
Replies
9
Views
2K
Top