Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Long Numerical Integration: Is this possible with MATLAB?

  1. Oct 4, 2011 #1
    Please, I am new in MATLAB and need some help as to whether MATLAB can perform this integration, and how to go about it. I have tried 'quad', 'quad8', 'trapz' and even a sum approach but it returns either a not-a-number answer or mtimes error or mrdivide error message. What I want to do is this (Pressure response of a horizontal observation well at a distance xd, yd and zd from a producer)
    integral= (1/(sqrt(td)))*(exp(yd^2/(4*td)))*(erf((1+xd)/(2*sqrt(td)))+erf((1-xd)/(2*sqrt(td))))*(1+2*(exp((-n^2*pi^2*td)/hd^2)*cos(n*pi*zwd)*cos(n*pi*((zd/hd)+zwd))));
    xd, yd, zd, hd and zwd are constants.
    i am integrating from 0 to td for td values of 0.001 to 10000
    for each of the integral segment, 0.001 to 10000, (exp((-n^2*pi^2*td)/hd^2)*cos(n*pi*zwd)*cos(n*pi*((zd/hd)+zwd))) will be summed for n=1 to infinity. ie Ʃ(exp((-n^2*pi^2*td)/hd^2)*cos(n*pi*zwd)*cos(n*pi*((zd/hd)+zwd))) for n=1 to ∞

    This is my recent attempt: where I did td =1 to 10000; and n= 1 to 100 just for it to work first.
    l=1000; %well half length
    % t=150;
    % td=(0.0002637*ky*t)/(phi*mu*ct*xf.^2)
    xd = x/l*sqrt(k/kx);
    yd = y/l*sqrt(k/ky);
    zd = z/h;
    hd = h/l*sqrt(k/kz);
    for i = 1:1e+1:10000
    for j=1:100
    int(i,j)= (1./(sqrt(td(i)))).*(exp(yd.^2./(4.*td(i)))).*(erf((1+xd)./(2.*sqrt(td(i))))+erf((1-xd)./(2.*sqrt(td(i))))).*(1+2.*(exp((-n(j)^2.*pi^2.*td(i))./hd.^2).*cos(n(j)*pi*zwd).*cos(n(j)*pi*((zd./hd)+zwd))));

    Error: ??? Undefined function or method 'td' for input arguments of type 'double'.
    **The attached document contains the equation.
  2. jcsd
  3. Oct 5, 2011 #2
    I'm having trouble understanding this. Is this all of your code? You don't have a definition for td.
  4. Oct 5, 2011 #3
    Yes, I am integrating with respect to td. So, I did not define td.
    Lower integral limit is zero, and I will have several upper limits of 0.001 to 10000 with varying intervals. The equation contains a sum series within.
    Equation is:
    integral= (1/(sqrt(td)))*(exp(yd^2/(4*td)))*(erf((1+xd)/(2*sqrt(td)))+erf((1-xd)/(2*sqrt(td))))*(1+2*Ʃ(exp((-n^2*pi^2*td)/hd^2)*cos(n*pi*zwd)*cos(n*pi*((zd/hd)+zwd))))
    the sigma is for n=1 to infinity.
  5. Oct 5, 2011 #4
    If you're using it as a symbolic variable, you need to define it using
    Code (Text):

    syms td;
    I don't know if there are alternate syntaxes for integrating, but AFAIK it's normally

    Code (Text):

    so your code would be

    Code (Text):

    syms td
    or something like that. If you're trying to integrate numerically rather than symbolically, then you should evaluate the function first and then integrate it, and you would need to define td.
  6. Oct 5, 2011 #5
    Thanks, let me try it.
Share this great discussion with others via Reddit, Google+, Twitter, or Facebook