Solving Fredholm Equation: Step-by-Step Guide

[huyhohoang]

Hi! I am now building a program to solve integral equation. I have understood the numerical method to solve it as well as others subroutines gauleg, ludcmp, lubksb but when I read this subroutine fred2 below, I don't know exactly the purpose of the author when putting the array omk(i,j) in the code. Can anyone give me some advice or instruction please?
Many thanks.

      SUBROUTINE fred2(n,a,b,t,f,w,g,ak)
      INTEGER n,NMAX
      REAL a,b,f(n),t(n),w(n),g,ak
      EXTERNAL ak,g
      PARAMETER (NMAX=200)
CU    USES ak,g,gauleg,lubksb,ludcmp
      INTEGER i,j,indx(NMAX)
      REAL d,omk(NMAX,NMAX)
      if(n.gt.NMAX) pause 'increase NMAX in fred2'
      call gauleg(a,b,t,w,n)
      do 12 i=1,n
        do 11 j=1,n
          if(i.eq.j)then
            omk(i,j)=1.
          else
            omk(i,j)=0.
          endif
          omk(i,j)=omk(i,j)-ak(t(i),t(j))*w(j)
11      continue
        f(i)=g(t(i))
12    continue
      call ludcmp(omk,n,NMAX,indx,d)
      call lubksb(omk,n,NMAX,indx,f)
      return
      END

 

[jedishrfu]

I'm not sure what the code is supposed to do but the omk array is setup basically like an identity matrix

and the line omk(i,j)-ak(t(i),t(j))*w(j) looks like the matrix equation OMK = I - AK*W

where I and AK are matrices and W is a vector.

OMK starts out in the code as the identity matrix and then its reused to hold the result of I - AK*W

so I guess the author was trying to save space by repurposing the OMK matrix memory.

 

[huyhohoang]

Thanks. Now I understand the purpose of the authors.