How can I solve an integral equation with an unknown kernel?

The equation I am trying to solve is:

where f(x) and h(x) are both complex and known, and g(x) is an unknown function. Presumably, the result should be a function g(x), however, it is not to be excluded that g(x) could actually be an operator instead. Can this be solved for either cases in MATLAB?

Thanks!

Answer 1

[Torsten]

g is not unique - it can be of any function type you like (we already had this discussion).

g(x)=1/integral_{x=0}^{x=2*pi} f(x)*h(x)dx

or

g(x)=1/(f(x)*h(x)*2*pi)

or

...

Best wishes

Torsten.

[Sergio Manzetti]

Thanks Torsten, :

with:

if 
y = cos(x).*e^(i.*x)
z = sin(x).*e^(i.*x)+sin(x).*cos(x)
G=1/integral_{x=0}^{x=2*pi} y.*z.*dx

>> Error: Line: 10 Column: 16 The expression to the left of the equals sign is not a valid target for an assignment.

[Torsten]

y=@(x)cos(x).*exp(1i*x);
z=@(x)sin(x).*exp(1i*x)+sin(x).*cos(x);
prod=y(x).*z(x);
G=1/integral(prod,0,2*pi);

Best wishes

Torsten.

[Sergio Manzetti]

There is something odd here: I still get an error of that integral form: Error using integral (line 82) First input argument must be a function handle.

Error in Integral_eqn (line 8) G=1/integral(prod,0,2*pi);

[Torsten]

prod=@(x)y(x).*z(x);

Best wishes

Torsten.

[Sergio Manzetti]

Thanks Torsten, I tried your suggestion on another set of functions, g and h and got: 398.942280401433e-003

Warning: Reached the limit on the maximum number of intervals in use. Approximate bound on error is 3.0e-01. The integral may not exist, or it may be difficult to approximate numerically to the requested accuracy. > In integralCalc/iterateScalarValued (line 372) In integralCalc/vadapt (line 132) In integralCalc (line 75) In integral (line 88) In Integral_eqn (line 8) >>

[Torsten]

If f and/or h have singularities in [0:2*pi], you will have to account for them in a suitable definition of g.

Best wishes

Torsten.

[Sergio Manzetti]

Dear Torsten, I Have to revisit this topic, because I have not been able to find C with the given command:

if true
  % code
end
syms h g x C
h = 1
g = 5
y=@(x)(C - (exp(2.*g.*1i.*x./h));
z=@(x)(C - (exp(-2.*g.*1i.*x./h));
prod=@(x)y(x).*z(x);
W= integral(prod,0,2*pi)==1;

can you advise on what to add to this command to find C?

Thanks!

[Sergio Manzetti]

PS I would expect the answer C to contain the interval value in it, so be in the form of 1/pi or similar. Alternatively the command should read:

if true
  % code
end
syms h g x C
h = 1
g = 5
y=@(x)(C - (exp(2.*g.*1i.*x./h));
z=@(x)(C - (exp(-2.*g.*1i.*x./h));
prod=@(x)y(x).*z(x)*(1+x^2);
W= integral(prod,0,L)==1;

where L is considered as an arbitrary constant, and the result be in the form of C = 1/2L or similar. Is this possible with this command? I get an error.

[Torsten]

Why don't you start from a solution that worked ?

syms L C x 
assume (L>0);
h = 1;
g = 5;
y = C-exp(2*g*1i*x/h);
z = C-exp(-2*g*1i*x/h);
prod = y*z*(1+x^2);
Csol = solve(int(prod,x,0,L)-1==0,C);

Best wishes

Torsten.

[Sergio Manzetti]

Dear Torsten, I did! However, you used a different solution, with "assume" in it I will try it now!

Thanks!

Answer 2

[John D'Errico]

If g(x) is unknown, then if all you have is a single equation equal to a constant, then there is no simple solution. Or, you can look at it as if there are infinitely many solutions, one of them being a constant function.

Just compute the integral of h(x)*f(x). Take the reciprocal. That is the value of the constant g that will make int(h*f*g) equal 1. So as long as int(h*f) over [0,2*pi] is not identically 0, then A solution is trivial. Yes there may be infinitely many other solutions, but they cannot be found unless you have information as to the functional form of g(x).

[Sergio Manzetti]

Thanks, I think the solution is not trivial, because the integral of h and f is:

W =

(150161193738355954586222597234585017169743443289080649840382803820382539880280999645654803891496677831582301321288235091356610775513114160448451324201927488546281633057*pi)/235872651551346334515097161311766805682307576276027526773512893765351510149509303255382920323603672971341773439422740306235661539370299539535089640714543502728087207936 - (857667986409676946755091410422420433765979555592970194729842043977566558535980733470902042660694127536712947665701541348018841649893734864216952836677423460461115319701245*sin((140170517671301*pi)/616503364953199))/27309427497697258096319858537523969165745773917815361139146936050648948950585748747769209155925329298757330388892969600751474273592254457836547309903960103405916641989001029353472 + (2731497041226781970742487761996229018221209713*pi*cos((140170517671301*pi)/616503364953199))/615656346818663737691860001564743965704370926101022604186692084441339402679643915803347910232576806887603562348544 + (13496744475957601245584225453924415705057921\\\n334313386459527769*pi^3)/6373655913801205184275058478629623527955117412724383078393329658404842718330226912869950027182751428026671953834958012184158952051832486872914239745199048223325975250024718628380737536

However, if I use "double(W) it is simply:

ans =

    1.99999997942694e+000

and the reciprocal is:

    500.000005143264e-003

and I have no information of the functional of g, except that it should be in the form e^(-i2*pi*x) and some constant. The absolute best way would be to solve g as an operator, and not as a function. That would solve this on a general level.

[Torsten]

Then make the Ansatz

g(x) = exp(-2*pi*i*x)+a

and determine a.

Best wishes

Torsten.

[Sergio Manzetti]

Yes, but John's suggestion was actually far easier, and it resulted that g has to be exactly 0.5 to make that integral go up to unity.

[John D'Errico]

Why is that non-simple? All it is is a constant. Just that the symbolic toolbox reports it as a mess. IT STILL HAS A RECIPROCAL!!!!!!!!! All numbers do, as long as they are not zero.

You are making this way more complex than it needs to be. It does not need to be a function. In fact, you have insufficient information to derive a functional form.

g = 1/W;

is entirely sufficient.

[Sergio Manzetti]

John: Did you read my comment at all? I wrote that your suggestion is more simple and it is only a constant, in agreement with your reply.

[Torsten]

But you wrote that g had to be of the form "e^(-i2*pi*x) and some constant". So you will have to use the above ansatz for g, don't you ?

Best wishes

Torsten.

[Sergio Manzetti]

That was a guess. But I suspect this in a non-trivial problem, because the functions (both of them) are discontinuous in R but continuous in C. So I have to go back to the drawing board and find what the original matrices are really "doing" in Hilbert space, OR I can make MATLAB numerically find an operator which makes these functions continuous in the space 0-2pi as a minimum space ), and make them equal to unity. I take that is not possible?

[John D'Errico]

Sorry. I missed the comment of yours that making g a constant did solve the problem. Since Answers only shows the last three comments by default, it is easy to miss one.

If you think about this as a linear functional on the family of functions g, so expand g as an orthogonal polynomial series in x, then one of those terms is a constant. And there is simply insufficient information to solve for more than one term. You could in theory solve for any function g that was known to within one unknown parameter. But a constant is sufficient.

[Sergio Manzetti]

Thanks John, is the "normalization" constant, or this constant (or a even a function) valid for these discontinuous systems treated without any modifications?

If we take for instance the Born-Sommerfeld condition, where the integral of two orthogonal functions is normalized to one, then the normalization constant is quite trivial to find. However, if the two functions are not orthogonal and not hermitian, they become discontinous and give complex values in the integral, when I use "double(W)" in MATLAB to determine their integral. Therefore, when if the discontinuous integral is normalized to one and g is expanded as a polynomial, would it work when h and f are not orthogonal?

[John D'Errico]

I don't see why not. Integration is just a linear operator. If g is a constant, then it can be pulled outside the integral.

[Sergio Manzetti]

Thanks I will give it another look.

Answer 3

[Sergio Manzetti]

Hi John, I found out that the non-hermitian nature of those functions makes them require a kernel, and not a constant, it appears to be (1+x^2). So I made this version to find the normalization constant when the integral is unity:

if true
  % code
end
syms h g x C
h = 1
g = 5
y=@(x)(C - (exp(2.*g.*1i.*x./h));
z=@(x)(C - (exp(-2.*g.*1i.*x./h));
prod=@(x)y(x).*z(x)*(1+x^2);
W= integral(prod,0,2*pi)==1;

and I should get the result for "C", however that is simply 0. It looks incorrect, or is the W part correctly written?

Thanks!