OU blog

Personal Blogs

Self ex-spline-atory

Visible to anyone in the world
Edited by Valentin Fadeev, Thursday, 27 Mar 2014, 10:06

M832 did not really fit in the "big picture" of my study plan this year, not only because I am notoriously bad at numerical calculations. Having invested a lot of effort and sleepless nights in developing intuition for the behaviour of analytic functions I was suddenly confronted by the cubic splines which seemed to have all those properties, that well-mannered functions would be never allowed to possess.

Nicely, but artificially glued together of several pieces of cubics, smooth only up to the second derivative, vanishing on the entire intervals, now this is what seems really counter-intuitive. Be that as it may, the TMA deadlines had to be met.

The following is the kind of problem I got stuck with for a while. Obviously I am not replicating a TMA question here, giving an extended solution to one of the problems from the Course Notes. So the task is to express a function, say f of x equals x squared in terms of cubic B-splines on the entire real axis. I am omitting a lot of background material focusing on one particular idea that arises in the solution.

equation left hand side x squared equals right hand side n ary summation from p equals negative normal infinity to normal infinity over lamda sub p times cap b sub p of x

Since cap b sub p cubed of x has a supporting interval equation left hand side open p comma sum with, 3 , summands p plus three plus one close equals right hand side open p comma p plus four close of length 4 outside which it vanishes, we can start by expressing the function in terms of cap b sub negative three comma cap b sub negative two comma cap b sub negative one comma cap b sub zero on open 0.1 close and then try to extend the result. Calculating the expressions for the splines on open zero comma one close :

cap b sub negative three of x equals negative one divided by 24 times open x minus one close cubed

cap b sub negative two of x equals one divided by 24 times open three times x cubed minus six times x squared plus four close

cap b sub negative one of x equals one divided by 24 times open sum with, 4 , summands negative three times x cubed plus three times x squared plus three times x plus one close

cap b sub zero of x equals one divided by 24 times x cubed

multiplying by the respective coefficients, summing and equating powers of x on each side we arrive at the following system of equations:

case statement case 1column 1 plus plus minus minus plus plus minus minus lamda minus minus three times times three lamda minus minus two times times three lamda minus minus one lamda zero equals equals zero case 2column 1 plus plus minus minus times times three lamda minus minus three times times six lamda minus minus two times times three lamda minus minus one equals equals zero case 3column 1 plus plus minus minus times times three lamda minus minus three times times three lamda minus minus one equals equals 24 case 4column 1 plus plus plus lamda minus minus three times times four lamda minus minus two lamda minus minus one equals equals zero

Having the solution (guaranteed by the Schoenberg-Whitney theorem): equation left hand side open lamda sub negative three comma lamda sub negative two comma lamda sub negative one comma lamda sub zero close equals right hand side open eight divided by three comma negative four divided by three comma eight divided by three comma 44 divided by three close

Now we want to find all coefficients on each of the intervals open xi sub p comma xi sub p plus four close for the points equation left hand side xi sub i equals right hand side i times h semicolon i equals prefix plus minus of one comma prefix plus minus of two full stop full stop full stop full stop . From the general expression for the B-spline it can be deduced that

cap b sub p of xi sub p plus one equals one divided by 24 times h

cap b sub p of xi sub p plus two equals one divided by six times h

cap b sub p of xi sub p plus three equals one divided by 24 times h

which for h equals one leads to the following recurrence relation:

equation left hand side sum with, 3 , summands lamda sub j minus one plus four times lamda sub j minus two plus lamda sub j minus three equals right hand side 24 times j squared

or, after changing the index

sum with, 3 , summands lamda sub j plus four times lamda sub j plus one plus lamda sub j plus two equals 24 left parenthesis j plus three times right parenthesis squared

Now here is the trick that I came up with and that was not (at least not explicitly) described in the Course Notes or the set book. The last expression can be thought of as a "second-order linear inhomogeneous recurrence relation". The advantage of this approach is that the structure of the solution instantly becomes clear.

The general solution of the corresponding homogeneous relation

sum with, 3 , summands lamda sub j plus four times lamda sub j plus one plus lamda sub j plus two equals zero

is derived in the course notes, using the standard method of solving this type of recurrencies and is given by the following expression:

equation left hand side lamda sub j super h equals right hand side alpha times open negative two minus Square root of three close super j plus beta times open negative two plus Square root of three close super j

It can also be found using generating functions. Not surprisingly it depends on 2 arbitary constants, as it takes 2 initial terms, lamda sub zero and lamda sub negative one to reconstruct the whole sequence from the three-term recurrency. Applying the general ideas from the linear systems we deduce that in order to obtain the general solution of the inhomogeneous recurrency we have to add a particular solution to the expression above.

Since the RHS is the quadratic polynomial it makes sence to look for the particular solution in the form:

equation left hand side lamda sub j super p equals right hand side sum with, 3 , summands a times j squared plus b times j plus c

Substituting this into the original recurrency and gatherig together the powers of j we obtain:

equation left hand side sum with, 5 , summands six times a times j squared plus open 12 times a plus six times b close times j plus eight times a plus six times b plus six times c equals right hand side sum with, 3 , summands 24 times j squared plus 144 times j plus 216

which after equating powers gives the solution equation left hand side open a comma b comma c close equals right hand side open four comma 16 comma 44 divided by three close

Thus the general solution of the inhomogeneous equation is given by the following formula:

equation left hand side lamda sub j equals right hand side sum with, 5 , summands alpha times open negative two minus Square root of three close super j plus beta times open negative two plus Square root of three close super j plus four times j squared plus 16 times j plus 44 divided by three

Now we can use the values of lamda sub zero and lamda sub negative one to determine the constants (bearing in mind that equation left hand side open negative two minus Square root of three close times open negative two plus Square root of three close equals right hand side negative one ):

equation left hand side 44 divided by three equals right hand side sum with, 3 , summands alpha plus beta plus 44 divided by three

equation left hand side eight divided by three equals right hand side negative alpha times open negative two plus Square root of three close minus beta times open negative two minus Square root of three close plus eight divided by three

which gives equation sequence alpha equals beta equals zero . Thus finally:

equation left hand side lamda sub j equals right hand side sum with, 3 , summands four times j squared plus 16 times j plus 44 divided by three

which is the solution of the original problem.

Permalink Add your comment
Share post

Around the infinity

Visible to anyone in the world
Edited by Valentin Fadeev, Sunday, 18 Sep 2011, 23:23

Now I got really fascinated with this topic. The most exciting part of evaluating integrals using residues is constructing the right contour. There are general guidelines for certain types of integrals, but in most cases the contour has to be tailored for a particular problem. Here is another example. Evaluate the following integral:

cap j equals integral over z sub one under z sub two Square root of sum with, 3 , summands cap a plus two times cap b divided by z plus cap c divided by z squared d z

Where cap a comma cap b comma cap c comma z sub one comma z sub two are real, cap a less than zero comma zero less than z sub one less than z sub two semicolon z sub one comma z sub two
are roots of the integrand.

The approach used in the previous example does not work as the integral along a large circle centered at the origin does not tend to 0. This, in fact, is a clue to the solution, as it prompts to have a look what indeed is happening to the function for large z . But first consider the integrand in the neighborhood of the origin (simple pole):

equation sequence Square root of sum with, 3 , summands cap a plus two times cap b divided by z plus cap c divided by z squared equals one divided by z times Square root of cap a times open z sub one minus z close times open z sub two minus z close equals one divided by z times Square root of cap a times z sub one times z sub two times open one minus z divided by z sub one close super one divided by two times open one minus z divided by z sub two close super one divided by two equals
equation sequence equals one divided by z times Square root of cap a times z sub one times z sub two times open one minus z divided by two times z sub one plus cap o of z squared close times open one minus z divided by two times z sub two plus cap o of z squared close equals one divided by z times Square root of cap a times z sub one times z sub two times open one minus z sub one plus z sub two divided by two times z sub one times z sub two times z plus cap o of z squared close
equation sequence equals one divided by z times Square root of cap a times cap c divided by cap a times open sum with, 3 , summands one plus cap b times cap a divided by cap a times cap c times z plus cap o of z squared close equals one divided by z times Square root of cap c times open sum with, 3 , summands one plus cap b divided by cap c times z plus cap o of z squared close
cap r sub zero of Square root of sum with, 3 , summands cap a plus two times cap b divided by z plus cap c divided by z squared equals Square root of cap c

Now the integrand is regular for large absolute value of z hence it can be expanded in the Laurent series convergent for cap r less than absolute value of z less than normal infinity where cap r is large:
equation sequence Square root of sum with, 3 , summands cap a plus two times cap b divided by z plus cap c divided by z squared equals Square root of cap a times open z sub one minus z close times open z sub two minus z close divided by z squared equals Square root of cap a times open one minus z sub one divided by z close super one divided by two times open one minus z sub two divided by z close super one divided by two equals

equation sequence equals Square root of cap a times open one minus z sub one divided by two times z plus cap o of one divided by z squared close times open one minus z sub two divided by two times z plus cap o of one divided by z squared close equals Square root of cap a times open one minus z sub one plus z sub two divided by two times one divided by z plus cap o of one divided by z squared close equals
equation left hand side equals right hand side Square root of cap a times open sum with, 3 , summands one plus cap b divided by cap a times one divided by z plus cap o of one divided by z squared close

The residue at infinity is defined to be the coefficient at one divided by z with the opposite sign:

cap r sub normal infinity of Square root of sum with, 3 , summands cap a plus two times cap b divided by z plus cap c divided by z squared equals negative cap b divided by Square root of cap a

Now we construct contour C by cutting the real axis along the segment open z sub one comma z sub two close and integrating along the upper edge of the cut in the negative direction, then along a small circle equation left hand side cap c sub z sub one equals right hand side open z colon absolute value of z minus z sub one equals epsilon close and along the lower edge of the cut. As the result of that one of the factors of the integrand increases its argument by equation left hand side e super two times pi times i divided by two equals right hand side e super pi times i
. Finally integrate along a small circle equation left hand side cap c sub z sub two equals right hand side open z colon absolute value of z minus z sub two equals epsilon close

cap j equals sum with, 4 , summands integral over z sub two minus epsilon under z sub one plus epsilon f of z d z plus integral over cap c sub z sub one f of z d z plus e super pi times i times integral over z sub one plus epsilon under z sub two minus epsilon f of z d z plus integral over cap c sub z sub two f of z d z
equation sequence lim over epsilon right arrow zero of integral over cap c sub z sub one f of z d z equals zero times lim over epsilon right arrow zero of integral over cap c sub z sub two f of z d z equals zero
lim over epsilon right arrow zero of cap j equals negative two times integral over z sub one under z sub two f of z d z

Integrating along this contour amounts to integrating in the positive direction along a very large circle centered at infinity. Hence the outside of C is the inside of this large circle which therefore includes both the residue at the origin and at infinity. Therefore, by the residue theorem:

cap j equals two times pi times i times open cap r sub zero plus cap r sub normal infinity close

equation sequence integral over z sub one under z sub two Square root of sum with, 3 , summands cap a plus two times cap b divided by z plus cap c divided by z squared d z equals negative pi times i times open Square root of cap c minus cap b divided by Square root of cap a close equals pi times i times open cap b divided by Square root of cap a minus Square root of cap c close

Permalink Add your comment
Share post

How comple^x can you get? Continued

Visible to anyone in the world
Edited by Valentin Fadeev, Thursday, 27 Mar 2014, 10:11

As a follow up thought I realized just how much easier it would have been to calculate the residue by definition, i.e. expanding the integrand in the Laurent series to get the coefficient a sub negative one . Let x equals negative one plus xi where xi is small:

equation sequence root of order four over four divided by left parenthesis one plus x times right parenthesis cubed equals root of order four over four divided by xi cubed equals two times e super pi times i divided by four divided by root of order four over four times left parenthesis one minus xi times right parenthesis super one divided by four left parenthesis one minus xi divided by two times right parenthesis super three divided by four divided by xi cubed

Of course, care is needed when choosing the value of the root. It depends on the value of the argument set on the upper bound of the cut. Since I chose it to be 0, the correct value of the root is e super pi times i divided by four .

Now

left parenthesis equation left hand side one minus xi times right parenthesis super one divided by four equals right hand side one minus one divided by four times xi minus three divided by 32 times xi squared plus cap o of xi cubed

left parenthesis equation left hand side one minus xi divided by two times right parenthesis super three divided by four equals right hand side one minus three divided by eight times xi minus three divided by 128 times xi squared plus cap o of xi cubed

left parenthesis one minus xi times right parenthesis super one divided by four left parenthesis equation left hand side one minus xi divided by two times right parenthesis super three divided by four equals right hand side one minus five divided by eight times xi minus three divided by 128 times xi squared plus cap o of xi cubed

Therefore, near x equals negative one the integrand has the following expansion:

equation left hand side root of order four over four divided by left parenthesis one plus x times right parenthesis cubed equals right hand side two times e super pi times i divided by four divided by root of order four over four times open one divided by left parenthesis x plus one times right parenthesis cubed minus five divided by eight left parenthesis x plus one times right parenthesis squared minus three divided by 128 times open x plus one close plus g of x close

where g of x is the regular part of the expansion which is of no interest in this problem.

Hence by definition:

equation sequence cap r sub negative one times f of x equals two times e super pi times i divided by four divided by root of order four over four times open negative three divided by 128 close equals negative three times e super pi times i divided by four divided by 64 times root of order four over four

Permalink
Share post

How comple^x can you get?

Visible to anyone in the world
Edited by Valentin Fadeev, Thursday, 21 Apr 2011, 23:58

I found this example in a textbook dated 1937 which I use as supplementary material for M828. It gave me some hard time but finally sorted out many fine tricks of contour integration. Some inspiration was provided by the discussion of the Pochhammer's extension of the Eulerian integral of the first kind by Whittaker and Watson.

Evaluate the following integral:

integral over zero under one root of order four over four divided by open one plus x close cubed d x

Consider the integral in the complex plane:

cap j equals contour integral over cap c root of order four over four divided by open one plus z close cubed d z

where contour C is constructed as follows. Make a cut along the segment open zero comma one close of the real axis. Let arg of z equals zero on the upper edge of the cut. Integrate along the upper edge of the cut in the positive direction. Follow around the point z equals one along a small semi-circle equation left hand side cap c sub gamma super one equals right hand side open z colon absolute value of one minus z equals gamma comma normal black letter cap i times z greater than zero close in the clockwise direction. The argument of the second factor in the numerator will decrease by three times pi divided by four . Then proceed along the real axis till z equals cap r and further round the circle equation left hand side cap c sub cap r equals right hand side open z colon absolute value of z equals cap r close in the counter-clockwise direction.

This circle will enclose both branch point of the integrand z equals zero and z equals one , however since the exponents add up to unity the function will return to its initial value:

equation sequence z super one divided by four times e super two times pi times i times one divided by four times open one minus z close super three divided by four times e super two times pi times i times three divided by four equals z super one divided by four times open one minus z close super three divided by four times e super two times pi times i equals z super one divided by four times open one minus z close super three divided by four

This is the reason why we only need to make the cut along the segment open zero comma one close and not along the entire real positive semi-axis. Then integrate from z equals cap r in the negative direction, then along the second small semi-circle equation left hand side cap c sub gamma squared equals right hand side open z colon absolute value of one minus z equals gamma comma normal black letter cap i times z less than zero close around the point z equals one where the argument of the second factor will again decrease by three times pi divided by four . Finally integrate along the lower edge of the cut and along a small circle equation left hand side cap c sub delta equals right hand side open z colon absolute value of z equals delta close around the origin where the argument of the first factor will decrease by equation left hand side two times pi times one divided by four equals right hand side pi divided by two .

contour

As the result of this construction the integral is split into the following parts:

multiline equation row 1 cap j equals sum with, 4 , summands integral over delta under one minus gamma f of z d z plus integral over cap c sub gamma super one f of z d z plus e super negative three times pi times i divided by four times integral over one plus gamma under cap r f of z d z plus integral over cap c sub cap r f of z d z plus row 2 sum with, 4 , summands prefix plus of e super five times pi times i divided by four times integral over cap r under one plus gamma f of z d z plus integral over cap c sub gamma squared f of z d z plus e super pi times i divided by two times integral over one minus gamma under delta f of z d z plus integral over cap c sub delta f of z d z

Two integrals around the small circles add up to an integral over a circle:

equation sequence integral over cap c sub gamma super one f of z d z plus integral over cap c sub gamma squared f of z d z equals integral over cap c sub gamma f of z d z equals integral over zero under two times pi root of order four over four divided by open two minus gamma times e super i times theta close cubed times i times gamma times e super i times theta d theta

lim over gamma right arrow zero of integral over cap c sub gamma f of z d z equals zero

Similarly, the integral around a small circle around the origin vanishes:

equation left hand side integral over cap c sub delta f of z d z equals right hand side integral over zero under two times pi root of order four over four divided by open one plus delta times e super i times theta close cubed times i times delta times e super i times theta d theta

lim over delta right arrow zero of integral over cap c sub delta f of z d z equals zero

Integrals along the segment open one plus gamma comma cap r close cancel out:

multiline equation row 1 equation left hand side e super negative three times pi times i divided by four times integral over one plus gamma under cap r f of z d z plus e super five times pi times i divided by four times integral over cap r under one plus gamma f of z d z equals right hand side e super negative three times pi times i divided by four times integral over cap r under one plus gamma f of z d z plus e super two times pi minus three times pi times i divided by four times integral over cap r under one plus gamma f of z d z equals row 2 equation sequence equals e super negative three times pi times i divided by four times integral over cap r under one plus gamma f of z d z minus e super negative three times pi times i divided by four times integral over cap r under one plus gamma f of z d z equals zero

The integral over the large circle also tends to 0 as R increases. This can be shown using the Jordan's lemma, or by direct calculation:

equation sequence integral over cap c sub cap r f of z d z equals integral over zero under two times pi root of order four over four divided by open one plus cap r times e super theta close cubed times i times cap r times e super i times theta d theta equals integral over zero under two times pi root of order four over four divided by open one divided by cap r plus e super i times theta close cubed times i times e super i times theta divided by cap r d theta

lim over cap r right arrow normal infinity of integral over cap c sub cap r f of z d z equals zero

Finally the only two terms that survive allow us to express the contour integral in terms of the integral along the segment of the real axis:

multiline equation row 1 equation sequence cap j equals integral over zero under one f of z d z minus e super pi times i divided by two times integral over zero under one f of z d z equals open one minus e super pi times i divided by two close times integral over zero under one f of z d z equals e super pi times i divided by four times open e super negative pi times i divided by four minus e super pi times i divided by four close times integral over zero under one f of z d z equals row 2 equation left hand side equals right hand side negative two times i times e super pi times i divided by four times sine of pi divided by four times integral over zero under one f of z d z

The contour cap c encloses the only singularity of the integrand which is the pole of the third order at z equals negative one . Hence, by the residue theorem:

cap j equals two times pi times i times cap r sub negative one times f of z

equation left hand side integral over zero under one f of z d z equals right hand side negative Square root of two times pi times e super negative pi times i divided by four times cap r sub negative one times f of z

The residue can be calculated using the standard formula:

equation sequence cap r sub negative one times f of z equals one divided by two times lim over z right arrow negative one of open open one plus z close cubed times root of order four over four divided by open one plus z close cubed close super double prime equals one divided by two times lim over z right arrow negative one of open root of order four over four close super double prime

Calculation of the derivative can be facilitated by taking logarithm first:

g equals root of order four over four

natural log of g equals one divided by four times natural log of z plus three divided by four times natural log of one minus z

equation sequence g super prime divided by g equals one divided by four times z minus three divided by four times open one minus z close equals one minus four times z divided by four times z times open one minus z close

equation left hand side g super prime equals right hand side one minus four times z divided by four times z super three divided by four times open one minus z close super one divided by four

natural log of g super prime equals natural log of one divided by four plus natural log of one minus four times z minus three divided by four times natural log of z minus one divided by four times natural log of one minus z

equation sequence g super double prime divided by g super prime equals negative four divided by one minus four times z minus three divided by four times z plus one divided by four times open one minus z close equals negative three divided by four times open one minus four times z close times z times open one minus z close

equation left hand side g super double prime equals right hand side negative three divided by 16 times z super seven divided by four times open one minus z close super five divided by four

equation sequence cap r sub negative one times f of z equals one divided by two times lim over z right arrow negative one of negative three divided by 16 times z squared times open one minus z close times z super one divided by four divided by open one minus z close super one divided by four equals negative three times e super pi times i divided by four divided by 64 times root of order four over four

equation sequence integral over zero under one root of order four over four divided by open one plus x close cubed d x equals negative Square root of two times pi times e super negative pi times i divided by four of negative three times e super pi times i divided by four divided by 64 times root of order four over four equals three times pi times root of order four over four divided by 64

Permalink
Share post

Going beyond dx

Visible to anyone in the world
Edited by Valentin Fadeev, Sunday, 18 Sep 2011, 23:23

This is quite a minor trick and like many things listed here may seem quite trivial. However, this is one of those few occasions when I had the tool in mind, before I actually got the example touse it on. Consider:

equation left hand side d times y divided by d times x equals right hand side one plus two divided by x plus y

which does not really require a great effort to solve. But forget all the standard ways for a moment and add d times x divided by d times x equals one to both parts:

equation left hand side d times open x plus y close divided by d times x equals right hand side two times one plus open x plus y close divided by x plus y

equation left hand side sum with, 3 , summands negative one plus one plus open x plus y close times d times open x plus y close divided by one plus open x plus y close equals right hand side two times d times x

equation left hand side d times open sum with, 3 , summands one plus x plus y close divided by one plus open x plus y close equals right hand side d times open y minus x close

natural log of absolute value of sum with, 3 , summands one plus x plus y equals open y minus x close plus natural log of cap c

equation left hand side sum with, 3 , summands one plus x plus y equals right hand side cap c times exp of y minus x

Hope this can be stretched to use in more complicated cases

Permalink Add your comment
Share post

Vector horror

Visible to anyone in the world
Edited by Valentin Fadeev, Sunday, 23 Jan 2011, 21:49

Had to do some revision of vector calculus/analysis before embarking on M828.

One point which I was not really missing, but did not quite get to grips with was the double vector product. I remembered the formula:

equation left hand side a right arrow multiplication open b right arrow multiplication c right arrow close equals right hand side b right arrow times open a right arrow dot operator c right arrow close minus c right arrow times open a right arrow dot operator b right arrow close ,

but nevertheless had difficulties applying it in excercises.

The reason whas that that the proof I saw used the expression of vector product in coordinates and comparison of both sides of the equation. However, I was aware of another, purely "vector" argument with no reference to any coordinate system.

Eventually I was able to reproduce only part of it, consulting one old textbook for some special trick. So here's how it goes.

For b right arrow multiplication c right arrow is perpendicular to the plane of b right arrow and c right arrow , a right arrow multiplication open b right arrow multiplication c right arrow close must lie in this plane, therefore:

equation left hand side a right arrow multiplication open b right arrow multiplication c right arrow close equals right hand side lamda times b right arrow plus mu times c right arrow

Dot-multiply both parts by a right arrow :

equation left hand side open a right arrow multiplication open b right arrow multiplication c right arrow close dot operator a right arrow close equals right hand side lamda times open b right arrow dot operator a right arrow close plus mu times open c right arrow dot operator a right arrow close

Since a right arrow multiplication open b right arrow multiplication c right arrow close up tack a right arrow , left-hand side is 0, so:

lamda times open b right arrow dot operator a right arrow close plus mu times open c right arrow dot operator a right arrow close equals zero reverse solidus qquod open asterisk operator close

Now define vector c prime right arrow lying in the plane of b right arrow and c right arrow , perpendicular to c right arrow and directed so that c prime right arrow , c right arrow and b right arrow multiplication c right arrow form the left-hand oriented system. This guarantees that the angle between b right arrow and c prime right arrow , open times times b right arrow c prime right arrow hat close less than pi divided by two .

Dot-multiply both parts by c prime right arrow :

equation left hand side open a right arrow multiplication open b right arrow multiplication c right arrow close dot operator c prime right arrow close equals right hand side lamda times open b right arrow dot operator c prime right arrow close