Re: beginner’s question with Maxi ma
Richard Fateman <[email protected]>
| Newsgroups | gmane.comp.mathematics.maxima.general |
|---|---|
| Message-ID | <CADB8Zm6GszyMdD3L=Se6-Yq3_Y_BOpbtmDknhD5FJtA+sZ=6FQ@mail.gmail.com> |
There are two separable parts to your computation. 1. The syntax hacking 2. The computation If you can refrain from changing the Maxima syntax and just use the given names and meanings of built-in functions, you should be able to state the computation and test it. That is, just use sum, product, floor, etc. If that works, you can play with the syntax to match what might appear to be natural for the application. For a user, for publication, for fun... Doing them both together, and, for instance, as was pointed out, redefining "-" can produce very obscure bugs. Good luck. RJF On Fri, Aug 21, 2026 at 5:19 PM <[email protected]> wrote: > Hello, all, > > I’m a new subscriber to this mailing list and an absolute beginner with > Maxima. I’m trying to use Maxima to calculate a meridian arc length on the > WGS 84 ellipsoid, and I seem to be missing something essential in my > understanding of Maxima. I’ve put my attempt so far into a .mac file, which > is presented below. I’d appreciate receiving your input on where my > knowledge has fallen short. Here’s my code in its current state: > > /* > meridian_arc_length.mac > > Calculate meridian arc lengths to zeptometer (zm) precision, > but display them at attometer (am) precision. > */ > fpprec: 21$ > fpprintprec: fpprec - 3$ > /* > Define the ⌊ ⌋ matchfix operator as an alias for floor(). > */ > matchfix("⌊", "⌋")$ > ⌊n⌋ := floor(n)$ > /* > Define the ², ³, and ⁴ postfix operators as aliases for ^2, ^3, and ^4 > respectively, with the same precedence as ^. > */ > exponent_precedence: 140$ > postfix("²", exponent_precedence, any, any)$ > "²"(n) := n^2$ > postfix("³", exponent_precedence, any, any)$ > "³"(n) := n^3$ > postfix("⁴", exponent_precedence, any, any)$ > "⁴"(n) := n^4$ > /* > Define the °, ′, ″, ‴, and ⁗ postfix operators for converting > arcdegrees, > arcminutes, arcseconds, arcthirds, and arcfourths respectively to > radians. > The precedence of these operators is just below those of the * and / > operators, so that e.g. “1/2°” is “(1/2)°” rather than “1/(2°)”. > */ > arc_precedence: 110$ > π: bfloat(%pi)$ > postfix("°", arc_precedence, any, any)$ > "°"(ad) := ad * π/180$ > postfix("′", arc_precedence, any, any)$ > "′"(am) := am * π/(180*60)$ > postfix("″", arc_precedence, any, any)$ > "″"(as) := as * π/(180*60²)$ > postfix("‴", arc_precedence, any, any)$ > "‴"(aþ) := aþ * π/(180*60³)$ > postfix("⁗", arc_precedence, any, any)$ > "⁗"(af) := af * π/(180*60⁴)$ > /* > Define the ∏() and ∑() macros for product() and sum() respectively, > so that their expressions are specified last rather than first. > */ > declare("∏", alphabetic)$ > ∏(i, i_min, i_max, expr) ::= product(expr, i, i_min, i_max)$ > declare("∑", alphabetic)$ > ∑(i, i_min, i_max, expr) ::= sum(expr, i, i_min, i_max)$ > /* > Load the ezunits package so that the semi-major axis > is specifiable in meters. > */ > load("ezunits")$ > /* > Use the WGS 84 ellipsoid definitions for Earth’s semi-major axis > and third flattening. > */ > declare("-", alphabetic)$ > semi-major_axis: 6378137 ` m; > first_flattening: 1/bfloat(298.257223563)$ > third_flattening: first_flattening/(2 - first_flattening); > /* > Find the meridian arc length from the equator to 45° (North or South). > The meridian arc length should be a bit less than semi-major_axis/2 > on the WGS 84 ellipsoid. > */ > latitude: 45°; > /* > Adapt S(φ), which is Equation (15) in the paper > “A General Formula for Calculating Meridian Arc Length […]” > (https://www.gsi.go.jp/common/000062452.pdf), > to be a memoizing function in Maxima, > except with a, n, and j_maximum specified as additional parameters, > so that the function can be used with other ellipsoids also. > > QUESTION: Have I correctly adapted the paper’s Equation 15 > for Maxima? (Equation 15 has the equivalent of inf for j_maximum, > and I wanted to ensure that low j_maximum values were working > before trying inf with S[].) > */ > S[a, n, φ, j_maximum] := \ > a * (1 - n)² * (1 + n) \ > * ∑(j, 0, j_maximum, ∏(k, 1, j, (-n/(2*k) - n)²)) \ > * (φ + ∑(l, 1, 2*j, (sin(2*l*φ)/l)) \ > * ∏(m, 1, l, \ > (-n/(2*j + 2*(-1)^m*⌊m/2⌋) - n)^((-1)^m))); > for j_max: 0 thru 5 do S[semi-major_axis, third_flattening, latitude, > j_max]; > /* > QUESTION: For each call of S[] with j_max from 0 through 5 inclusive, > all of j through n are not completely resolved in the returned values. > > For example, with j_max = 0, > > (%i1) S[semi-major_axis, third_flattening, 45°, 0]; > 2 j > │-n │ > (%o1) (6.36740874757294245b6 ` m) │─── - 1.67922038638370455b-3│ > │2 k │ > ((2 j sin(1.57079632679489662b0 l) > m > -n -1 l > ((──────────────────── - 1.67922038638370455b-3) ) )/l > m m > 2 -1 floor(─) + 2 j > 2 > + 7.8539816339744831b-1) > (%i2) _ > > What have I failed to do to ensure the resolution of j through n? > */ > > > _______________________________________________ > Maxima-discuss mailing list > [email protected] > https://lists.sourceforge.net/lists/listinfo/maxima-discuss > _______________________________________________ Maxima-discuss mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/maxima-discuss