Mathematica 9 is now available
Services & Resources / Wolfram Forums / MathGroup Archive
-----

MathGroup Archive 2011

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Integrating a Piecewise function four times but different results

  • To: mathgroup at smc.vnet.net
  • Subject: [mg121791] Re: Integrating a Piecewise function four times but different results
  • From: Peter Pein <petsie at dordos.net>
  • Date: Sun, 2 Oct 2011 02:37:08 -0400 (EDT)
  • Delivered-to: l-mathgroup@mail-archive0.wolfram.com
  • References: <j66e9t$hch$1@smc.vnet.net>

Am 01.10.2011 09:08, schrieb Rob Y. H. Chai:
> Hello everyone,
> 
>  
> 
> I tried to integrate a piecewise function 4 times but got different results
> for two supposedly identical cases. The integration was done in the Mathematica
> version 8.0.1.0. I am interested in what I did wrong, so any comments will
> be appreciated. My code is below.
> 
>  
> 
> Sincerely
> 
>  
> 
> Rob Y. H. Chai
> 
> UC Davis
> 
>  
> 
>  
> 
> (* Introduction p1[x] as piecewise funtion *)
> 
> p1[x_] = Piecewise[{{2*wmax*x/L, x <= L/2}, {2*wmax*(1 - x/L), x > L/2}}];
> 
>  
> 
>  
> 
> (* p2[x] is a specialized case of p1[x] with wmax = 1 and L = 1 *)
> 
> p2[x_] = p1[x] /. wmax -> 1 /. L -> 1; 
> 
>  
> 
>  
> 
> (* now integrate p1[x] 4 times and introduce integration constants *)
> 
> soln1[x_] = Integrate[Integrate[Integrate[Integrate[p1[x], x], x], x], x] +
> a*x^3 + b*x^2 + c*x + d ;
> 
>  
> 
> (* solve for integration constants in the solution of p1[x] *)
> 
> eq1a = soln1[0] == 0;
> 
> eq1b = (D[soln1[x], {x, 2}] /. x -> 0) == 0;
> 
> eq1c = soln1[L] == 0;
> 
> eq1d = ( D[soln1[x], {x, 2}] /. x -> L) == 0;
> 
> soln1A[x_] = soln1[x] /. Simplify[Solve[{eq1a, eq1b, eq1c, eq1d}, {a, b, c,
> d}], L > 0][[1]];
> 
> (* check case 1 solution at x = L/2 *)
> 
> soln1A[L/2] /. wmax -> 1 /. L -> 1
> 
> 17/1920
> 
>  
> 
>  
> 
> (* similarly integrate p2[x] 4 times and introduce integration constants *)
> 
> soln2[x_] = Integrate[Integrate[Integrate[Integrate[p2[x], x], x], x], x] +
> r*x^3 + s*x^2 + u*x + v;
> 
> (* solve for integration constants in the solution of p2[x] *)
> 
> eq2a = soln2[0] == 0;
> 
> eq2b = (D[soln2[x], {x, 2}] /. x -> 0) == 0;
> 
> eq2c = soln2[1] == 0;
> 
> eq2d = ( D[soln2[x], {x, 2}] /. x -> 1) == 0;
> 
> soln2A[x_] = soln2[x] /. Simplify[Solve[{eq2a, eq2b, eq2c, eq2d}, {r, s, u,
> v}]][[1]];
> 
> soln2A[1/2] (* case 2 solution at x = 1/2 *)
> 
> 1/120
> 
>  
> 
> Plot[{soln1A[x] /. wmax -> 1 /. L -> 1, soln2A[x]}, {x, 0, 1}]
> 
>  
> 
> Question - why doesn't case 1 solution converge to case 2 solution when case
> 2 is simply a special case of case 1?

Case 1 has got a discontinuity at L/2 (already after the first
integration). I would prefer to iterate (Integrate w.r.t. x from 0 to z,
add a constant and replace z by x) four times:

assum = {L > 0 < wmax, Element[z, Reals]};

{p0[x_], p2[x_]} =
   FoldList[FullSimplify[PiecewiseExpand[
      c[#2] + Integrate[#1, {x, 0, z}, Assumptions -> assum] /. z -> x],
      assum] & ,
    Piecewise[{{(2*wmax*x)/L, x <= L/2}, {2*wmax*(1 - x/L), x > L/2}}],
    Range[4]]
 [[{-1, -3}]];

will give us the candidate function p0 and its second derivative.

the conditions are:

eqns = Simplify[Flatten[(#1 /@ {0, L} & ) /@ {p0, p2}] == 0, assum];

so we get:

In[4]:= pfunc[x_] = FullSimplify[p0[x] /. First[
          Solve[eqns, Union[Cases[eqns, _c, Infinity]]]],
        assum]

Out[4]= Piecewise[{
         {(wmax*x*(5*L^2 - 4*x^2)^2)/(960*L), x == 0 || L >= 2*x}},
        (wmax*(L - x)*(L^2 + 8*L*x - 4*x^2)^2)/(960*L)]

and the lim_{x->L/2}{pfunc(x)} is the same for both directions:

In[5]:= (Limit[pfunc[x], x -> L/2,
         Assumptions -> assum, Direction -> #1] & ) /@ {1, -1}

Out[5]= {(L^4*wmax)/120, (L^4*wmax)/120}

Hope that helps,
Peter



  • Prev by Date: A collection of Mathematica learning resources
  • Next by Date: Re: ParallelDo and C-compiled routines
  • Previous by thread: Re: Integrating a Piecewise function four times but different results
  • Next by thread: Re: Re: MakeExpression and color