Home |
Search |
Today's Posts |
|
#1
![]()
Posted to microsoft.public.excel.misc
|
|||
|
|||
![]() I've been using excel for some math analysis and came to a stop when I found that excel did not have an integration function I'm not asking for anything fancy. Numerical integration using the trapezoidal rule would do. I have some foggy ideas how to go about this. Here goes. Write the integration function in qbasic(a primitive language but still good enough) then tell qbasic to make a file containing the data that is needed. Then excel SOMEHOW looks at this file and pastes the data into the cells. Is there a better way?? or how does one implement the above method? -- integreat ------------------------------------------------------------------------ integreat's Profile: http://www.excelforum.com/member.php...o&userid=34282 View this thread: http://www.excelforum.com/showthread...hreadid=540492 |
#2
![]()
Posted to microsoft.public.excel.misc
|
|||
|
|||
![]()
I ran a search in this newgroup for integration.
Niek Otten provides a link in response to a question perhaps similar to yours. HTH http://www.microsoft.com/office/comm...195&sloc=en-us "integreat" wrote: I've been using excel for some math analysis and came to a stop when I found that excel did not have an integration function I'm not asking for anything fancy. Numerical integration using the trapezoidal rule would do. I have some foggy ideas how to go about this. Here goes. Write the integration function in qbasic(a primitive language but still good enough) then tell qbasic to make a file containing the data that is needed. Then excel SOMEHOW looks at this file and pastes the data into the cells. Is there a better way?? or how does one implement the above method? -- integreat ------------------------------------------------------------------------ integreat's Profile: http://www.excelforum.com/member.php...o&userid=34282 View this thread: http://www.excelforum.com/showthread...hreadid=540492 |
#3
![]()
Posted to microsoft.public.excel.misc
|
|||
|
|||
![]()
Hi. This should do it. This is something I haven't yet finished up. You
might also do a Find for "HOWDY" to go beyond trapezoidal. HTH Dave Braden A Worksheet-based Solution Suppose your x-values are in a column-range you named Xpts, and the corresponding y-values are in a columns-range Ypts. What is the area under the curve? To keep the following worksheet formula short, choose Insert/Name/Define, and set "Pnls" to be =COUNT(Xpts)-1. This is the number of panels formed by successive pairs of (x1,x2) points. Alternatively, you can enter the formula into a cell, and then Insert/Name/Define that cell to be Pnls. Next, enter this function into a convenient cell: =SUMPRODUCT((OFFSET(Xpts,1,0,Pnls)-OFFSET(Xpts,0,0,Pnls)),(OFFSET(Ypts,1,0,Pnls)+OFFS ET(Ypts,0,0,Pnls)))*0.5 The values in Xpts must be nondecreasing (e.g., 1.2, 2.3, 2.3, 3.0, 4.0) or all non-increasing (e.g., 5.2, 4.0, 4.0, 3.1, 2.5). You can build in some error-checking for this by modifying the above formula to =SUMPRODUCT((OFFSET(Xpts,1,0,Pnls)-OFFSET(Xpts,0,0,Pnls)),(OFFSET(Ypts,1,0,Pnls)+OFFS ET(Ypts,0,0,Pnls)))*0.5*IsWeakMonotone where "IsWeakMonotone" is the name of a cell containing the array-entered formula =IF(OR(AND(OFFSET(Xpts,1,0,Pnls)=OFFSET(Xpts,0,0, Pnls)),AND(OFFSET(Xpts,1,0,Pnls)<=OFFSET(Xpts,0,0, Pnls))),1,NA()) Remember, to enter an array formula, hold ctrl-shift when pressing Enter. This is a (theoretically) exact solution to the problem. It sums up the areas of the individual trapezoids formed by the panels. To get a graphic feel for this, and where we're headed below, try it! If Xpts forms a *row* range, adjust the corresponding OFFSETs above for Xpts to OFFSET(Xpts,0,1,0,Pnls) and OFFSET(Xpts,0,0,0,Pnls). Do similarly for Ypts as needed. If one is a row and the other a column, use Excel's TRANSPOSE function (see Excel's Help if you are unfamiliar with this); I'd recommend using it on the row-oriented range in the Named Definition. A VBA Solution If you want a more flexible solution, you can use the following VBA (Visual Basic for Applications) code. It checks for weak monotonicity of Xpts, and other potential problems. This code has been tested for Excel versions 9 and later; earlier versions of Excel may require some editing. To use, enter into a convenient cell =TrapezIntegrate(Xpts,Ypts) where each of the arguments is either a single-row or single-column range of cells that evaluate to numbers Optional arguments let you determine the area under part of the curve. Option Explicit '//////////////////////////////////////////////////////////////////////////////// ' TrapezIntegrate ' For given arrays of (x,y) pts, determine integral, ' assuming pts connected by straight lines, and that ' the x-values are nondecreasing. ' ' David J. Braden ' 2002 March 13 '//////////////////////////////////////////////////////////////////////////////// Function TrapezIntegrate(Xin As Variant, Yin As Variant, _ Optional LowX As Double = -1E+307, Optional HighX As Double = 1E+307 _ ) As Variant Dim lNumPts As Long, l As Long Dim lLoPart As Long, lHiPart As Long Dim bNegInt As Boolean Dim dDelta As Double, dSum As Double Dim Xpts, Ypts 'variants On Error Resume Next If Not MakeVect(Xin, Xpts) Then TrapezIntegrate = Xpts 'return w/ error Exit Function End If lNumPts = UBound(Xpts) If Not MakeVect(Yin, Ypts) Then TrapezIntegrate = Ypts 'return w/ error Exit Function ElseIf (lNumPts < UBound(Ypts)) Or (lNumPts < 2) Then TrapezIntegrate = CVErr(xlErrValue) Exit Function End If If LowX HighX Then 'switch them dDelta = HighX: HighX = LowX: LowX = dDelta: bNegInt = True End If If HighX < Xpts(1) Or LowX Xpts(lNumPts) Then TrapezIntegrate = CVErr(xlErrNum) Exit Function End If If LowX < Xpts(1) Then LowX = Xpts(1) lLoPart = 1 Else lLoPart = Application.WorksheetFunction.Match(LowX, Xpts, 1) End If If HighX Xpts(lNumPts) Then HighX = Xpts(lNumPts) lHiPart = lNumPts - 1 Else lHiPart = Application.WorksheetFunction.Match(HighX, Xpts, 1) If lHiPart = lNumPts Then lHiPart = lHiPart - 1 End If 'check that Xpta and Ypts numeric outside of range to be integrated 'the interior pots will be checked during integration With Application.WorksheetFunction For l = 1 To lLoPart If Not (.IsNumber(Xpts(l)) And .IsNumber(Ypts(l))) Then TrapezIntegrate = CVErr(xlErrValue) Exit Function End If Next For l = lHiPart + 1 To lNumPts If Not (.IsNumber(Xpts(l)) And .IsNumber(Ypts(l))) Then TrapezIntegrate = CVErr(xlErrValue) Exit Function End If Next End With If lLoPart = lHiPart Then 'points are in the same partition dSum = (LinTerp(LowX, Xpts(lLoPart), Xpts(lLoPart + 1), Ypts(lLoPart), Ypts(lLoPart + 1)) _ + LinTerp(HighX, Xpts(lLoPart), Xpts(lLoPart + 1), Ypts(lLoPart), Ypts(lLoPart + 1))) _ * (HighX - LowX) Else dSum = (LinTerp(LowX, Xpts(lLoPart), Xpts(lLoPart + 1), Ypts(lLoPart), Ypts(lLoPart + 1)) _ + Ypts(lLoPart + 1)) * (Xpts(lLoPart + 1) - LowX) _ + (LinTerp(HighX, Xpts(lHiPart), Xpts(lHiPart + 1), Ypts(lHiPart), Ypts(lHiPart + 1)) _ + Ypts(lHiPart)) * (HighX - Xpts(lHiPart)) For l = lLoPart + 1 To lHiPart - 1 dDelta = Xpts(l + 1) - Xpts(l) dSum = dSum + dDelta * (Ypts(l + 1) + Ypts(l)) If Err.Number < 0 Then 'one of the (X,Y) not numeric TrapezIntegrate = CVErr(xlErrValue) Exit Function ElseIf Sgn(dDelta) < 0 Then TrapezIntegrate = CVErr(xlErrNum) Exit Function End If Next End If If bNegInt Then TrapezIntegrate = CVar(-dSum * 0.5) Else TrapezIntegrate = CVar(dSum * 0.5) End If End Function '//////////////////////////////////////////////////////////////////////////////// Private Function LinTerp(x, x1, x2, y1, y2) On Error Resume Next If x1 = x2 Then If x = x1 Then LinTerp = y2 Else LinTerp = CVErr(xlErrNum) Else LinTerp = y1 + (y2 - y1) * (x - x1) / (x2 - x1) End If If Err.Number < 0 Then LinTerp = CVErr(xlErrValue) End Function '//////////////////////////////////////////////////////////////////////////////// ' Makevect takes an array or range, and orients it as a row vector ' Return error if vSrc has more than 1 row and more than 1 col. '//////////////////////////////////////////////////////////////////////////////// Public Function MakeVect(vSrc As Variant, vx As Variant) As Boolean Dim lCol As Long, lRow As Long If IsError(vSrc) Then vx = vSrc: Exit Function If TypeName(vSrc) = "Range" Then If vSrc.Areas.Count 1 Then vx = CVErr(xlErrValue): Exit Function vx = vSrc.Value Else vx = vSrc End If On Error Resume Next lRow = UBound(vx, 1): lCol = UBound(vx, 2) If lCol = 1 Then 'already have a column MakeVect = True vx = Application.WorksheetFunction.Transpose(vx) ElseIf lCol = 0 Then 'have a scalar or row MakeVect = True If lRow = 0 Then ReDim vx(1 To 1): vx(1) = vSrc Else 'have at least 2 columns If lRow 1 Then vx = CVErr(xlErrValue) Else With Application.WorksheetFunction vx = .Transpose(vx) vx = .Transpose(vx) End With MakeVect = True End If End If End Function B) The Area with Smoothing Microsoft Excel apparently uses Bezier smoothing for its charts. The downside of this, as far as the topic at hand goes, is that even if the support of the function (the "X points") is monotone (i.e., nondecreasing or nonincreasing), the curve might bend back on itself, rendering the question "what is the area under the curve" meaningless (or at best ambiguous). To see this, enter into A1:A5 the values 1, 2, 3, 3.1, 4.5, and into B1:B5 enter 1.5, 2, 1, 3, 2.5. Select A1:B5 and create a scatter (XY) chart. Notice how the curve "bends back" between 3 and 3.1. However, the Bezier smooth is but one of many. There are two in particular that warrant attention: the cubic spline, and the (double) parabolic. 1) Cubic Spline The cubic spline is very smooth; it is far better behaved than the Bezier in maintaining monotonicity, yet tends, by comparison, to "overshoot" under when, say, the data pairs take a sudden change, as in Xpts = {1,2,2.1,3}, Ypts = {2,2,5,3} 2. I Know the Function I am Plotting Easily the fastest way to do this within Excel is with worksheet functions, but see below for VBA solutions. Reference material follows. HOWDY A) Worksheet Solutions Here's a way to set up your worksheet to get three popular solutions, all of which fall in the Newton-Cotes family. They are the Trapezoidal Rule, Simpson's Rule, and Weddle's Rule. To use Excel for evaluating the integral of, say, 2+3*(Ln(x))^0.6: ¥ Enter some labels: In cell A1, enter "X_1" A2 "X_n" A3 "NbrPanels" ¥ Set some values: In cell B1, enter 1 B2 2.5 B3 600 ¥ Define some names: Select A1:B3, then choose Insert/Name/Create. Make sure that only the "Left Column" box is selected. If it isn't, you might have entered text instead of numbers in the right column, or mis-selected the range. Press OK. ¥ Choose Insert/Names/Define Enter each of the following names and their definitions, pressing Add with each entry (you can copy and paste these): EPanels =6*ROUNDUP(NbrPanels/6,0) delta =(X_n-X_1)/EPanels Steps =ROW(INDIRECT("1:"&EPanels+1))-1 EvalPts =X_1+delta*Steps WeddleWts =IF(MOD(Steps,EPanels)=0,1,IF(MOD(Steps,6)=0,2,IF( MOD(Steps,3)=0,6,IF(MOD(Steps,2)=0,1,5))))*(0.3*De lta) EPanels takes NbrPanels and rounds it up to the next multiple of 6. This isn't neccessary for Simpson's rule, which would need =EVEN(NbrPanels) instead, or the Trapezoidal, which needs no adjustment. It *is* required, though, for the more sophisticated Weddle's Rule. (optional) If you are interested in a trapezoidal approximation, define TrapWts =IF(MOD(Steps,EPanels)=0,0.5*Delta,Delta) For Simpsons's Rule, define SimpsonWts =IF(MOD(Steps,EPanels)=0,1,IF(MOD(Steps,2)=1,4,2)) *(Delta/3) ¥ Close the Define Names box, and in, say, cell D1, array-enter =SUM(WeddleWts*(2+3*LN(EvalPts)^0.6)) That is, type in the function, and hold Ctl-Shift when pressing Enter. In general, ctrl-shift-enter =SUM(MyWts*f(EvalPts)) where f() is a legitimate Excel expression that yields a scalar numeric value, and MyWts is any one of TrapWts, SimpsonWts or WeddleWts. For the sample problem, Mathematica gives, at a higher precision than Excel works with, 5.94249912704062. Notes: ¥ NbrPanels should be at least 6 for this implementation. You typically achieve greater accuracy as NbrPanels increases. ¥ You (usually) do better, for given NbrPanels, with SimpsonWts than TrapWts, and better yet with WeddleWts. I don't know of a function, *as graphed by Excel*, which isn't best handled by WeddleWts. ¥ If you have "jumps" in your function, break it up at the points where those occur, and add the pieces. ¥ To keep keep your worksheet so that the calculations update with changes in X_1, X_n or NbrPanels, enter =Delta in some cell (say, $B$4), and be sure that Calculation isn't set to Manula in Excel's Preferences. ¥ Using a high value for NbrPanels (greater than 6^4 = 1296 or so) might noticably slow down your workbook's performance under certain circumstances. Chances are, especially with WeddleWts, you will need far fewer points to get your desired accuracy. B) VBA Solutions While VBA solutions are perhaps slower than worksheet-based ones, they offer more flexibility. If you want a one-time solution, then the following is quick and easy. Paste the code below into a VBA module, and run CalcInt (note-this uses a slightly modified version of Karl's post): '////////////////// Option Explicit Sub CalcInt() Dim x As Double x = Romberg(0, 2.5, 13) Debug.Print " x= " & x End Sub ' Hard-wire your function here Function fct(x As Double) As Double fct = 2 + 3 * (Log(x) ^ 0.6) End Function ' The function Romberg calculates the value of the integral of ' the function fct() between the limits lgr and rgr. ' The order of the error term is 2*Order+2, Order=0. ' Karl M. Syring, 1/6/01 ) Function Romberg(LowX As Double, HiX As Double, Order As Long) As Double ReDim t(1 To Order + 1) As Double Dim dSup As Double, dU As Double, dM As Double Dim f As Long, h As Long, j As Long, n As Long dSup = HiX - LowX t(1) = (fct(LowX) + fct(HiX)) * 0.5 n = 2 For h = 1 To Order dU = 0# dM = dSup / n For j = 1 To n - 1 Step 2 dU = dU + fct(LowX + j * dM) Next j t(h + 1) = dU / n + t(h) * 0.5 f = 1 For j = h To 1 Step -1 f = 4 * f t(j) = t(j + 1) + (t(j + 1) - t(j)) / (f - 1) Next j n = 2 * n Next h Romberg = t(1) * dSup End Function '////////////////// Running the above (where the third argument of Romberg = 13) will give, for the particular example specified in "fct", a result comparable to the worksheet-based approaches using NbrPanels =2*6^5; the relative error is about 2E-08 among the methods. For more general situations (e.g., not needing to "hard-code" the function definition), and more advanced algorithms, contact the author of this note. 3. References Press, et al., Numerical Recipes in C, available on-line at http://lib-www.lanl.gov/numerical/bookcpdf.html Eric Weisstein's "World of Mathematics" web page integreat wrote: I've been using excel for some math analysis and came to a stop when I found that excel did not have an integration function I'm not asking for anything fancy. Numerical integration using the trapezoidal rule would do. I have some foggy ideas how to go about this. Here goes. Write the integration function in qbasic(a primitive language but still good enough) then tell qbasic to make a file containing the data that is needed. Then excel SOMEHOW looks at this file and pastes the data into the cells. Is there a better way?? or how does one implement the above method? -- Please keep response(s) solely within this thread. |
#4
![]()
Posted to microsoft.public.excel.misc
|
|||
|
|||
![]() I know how to do machine language programming so I guess it will take some patience to learn VBA Thanks for response -- integreat ------------------------------------------------------------------------ integreat's Profile: http://www.excelforum.com/member.php...o&userid=34282 View this thread: http://www.excelforum.com/showthread...hreadid=540492 |
#5
![]()
Posted to microsoft.public.excel.misc
|
|||
|
|||
![]()
You're welcome.
Don't miss, though, the worksheet-based solutions to the problem that I posted in that note. Dave B integreat wrote: I know how to do machine language programming so I guess it will take some patience to learn VBA Thanks for response -- Please keep response(s) solely within this thread. |
Reply |
Thread Tools | Search this Thread |
Display Modes | |
|
|
![]() |
||||
Thread | Forum | |||
Does True Have a Numerical Value | Excel Worksheet Functions | |||
Changing Y-Axis Labels from Numerical to Textual | Charts and Charting in Excel | |||
Return Numerical Label for LAST value Subtracted to reach Sum Target Value | Excel Worksheet Functions | |||
Numerical grade to Alpha character | Excel Discussion (Misc queries) | |||
Hiding Personal.xls using Novell Groupwise integration | Excel Discussion (Misc queries) |