I get this in Maple 2018.2 and Maple 2020.0 (64bit Linux). This is a shorter and slightly different example.

restart;
kernelopts(version);
    Maple 2020.0, X86 64 LINUX, Mar 4 2020, Build ID 1455132

sum1 := Sum(1/((11-2*n)!*(n-1)!*(29+2*n)!),n=2..3):

value(sum1) ;
                              173                    
          -------------------------------------------
          7439866535798024349359988963016704000000000

simplify(sum1);
Error, (in evalr/evalr) too many levels of recursion

foo:=convert(sum1,GAMMA);
              3                                           
            -----                                         
             \                                            
              )                      1                    
     foo :=  /    ----------------------------------------
            ----- GAMMA(12 - 2 n) GAMMA(n) GAMMA(30 + 2 n)
            n = 2                                         

combine(foo);
                               0

simplify(foo);
                               0

value(foo);
                              173                    
          -------------------------------------------
          7439866535798024349359988963016704000000000

simplify(sum1);
                               0

I'll point out that simplify(sum1) returned 0 at the end, although earlier it produced a recursion limit error. On some runs combine(foo) would also error out like that, with foo in GAMMA form.

simplify_ex_ac.mw

And one more short example (without the (n-1)! in the denominator). simplify_exac2.mw

Someone might mention that add is a better choice for adding up a finite number of terms. But that seems a mere sidebar here. The following is also interesting, where the case that N<5 is pushed aside.

restart; sum(1/((11-2*n)!*(29+2*n)!),n=2..N,formal=true);

                             3079                     
         ---------------------------------------------
         130088533681106143868879347831013376000000000

restart; sum(1/((11-2*n)!*(29+2*n)!),n=2..N,formal=false);

  Error, (in SumTools:-DefiniteSum:-ClosedForm) summand is
  singular in the interval of summation

And,

restart;
combine(Sum(1/(GAMMA(7-2*n)*GAMMA(1+2*n)),n=2..3,formal=false));
                               0

restart;
combine(Sum(1/(GAMMA(7-2*n)*GAMMA(1+2*n)),n=2..3,formal=true));
                               0

It looks like your plan is to have the module local atable be assigned the actual table, rather than just a name that refers to that table. For example,

TestLibrary := module()
    export aset,atable:
    option package:

    aset := myset:
    atable := eval(mytable);
end module

If -- as you had it -- the module local is only assigned the name of the table (ie. mytable) then that's what gets stored in the archive.

This evaluation behavior is special to a few data structures, including table and Record -- but not list, set, Matrix/Vector/Array.

You might want to read the Help pages for table , in particular the 5th bullet point which describes that a table has special evaluation rules. See also last name evaluation .

Attached is a revision that preserves the table, upon savelib and restart.
TestLibrary_ac.mw

You attached a data file, but you stil have not attached any worksheet showing your actual code.

In a few of your responses you included something like code, as text. But it was all incomplete, and with different sizes than the attached data, and very incomplete since none of it worked alone. And there were syntax errors in extracting particular columns.

Anyway, the following works for me. I used LinearAlgebra:-Dimensions to pick off the number of rows and columns of the imported data. I used square brackets to index into the Matrix -- over a range -- to extract a single column as Vector.

signalplot_ac.mw

Yes, there are problems here. For example, using Maple 2019.2 or Maple 2020.0,

restart;
assume(0 < a, a < R);

K := 2*Pi*Int(sin(phi)
              *Int(rho^2,
                   rho = sin(phi)*R - sqrt(sin(phi)^2*R^2 - R^2 + a^2)
                         ..
                         sin(phi)*R + sqrt(sin(phi)^2*R^2 - R^2 + a^2)),
              phi = arccos(a/R) .. Pi - arccos(a/R)):

value(K);

                -Pi^2*R*a^2

#
# The following call to `simplify` merely changes the radical
# in the inner limits of integration to,
#
#    ( -R^2*cos(phi)^2 + a^2 )^(1/2)
#
simplify(K, trig):
value(%);

                2*Pi^2*R*a^2

#
# The following call to `combine` merely changes the radical
# in the inner limits of integration to,
#
#    ( -2*R^2 + 4*a^2 - 2*R^2*cos(2*phi) )^(1/2)/2
#
combine(K):
value(%);

                2*Pi^2*R*a^2

Below we can seen that the inner integral is not the problem. The problem comes from the elliptictrig method which is one of those that int tries for the outer integral, and whose result -- when it does not fail -- gets preferred and returned.

One possible way to avoid the problem is to force int to use a method other than the problematic elliptictrig. This means not using the nicely typeset 2D Input integral, so I'll use plaintext 1D Maple Notation instead. Below I'll continue in the same session, with the same assumptions on a and R.

PP := int(rho^2,
          rho = sin(phi)*R - sqrt(-R^2*cos(phi)^2 + a^2)
                ..
                sin(phi)*R + sqrt(-R^2*cos(phi)^2 + a^2));

     1/3*(sin(phi)*R+(-R^2*cos(phi)^2+a^2)^(1/2))^3
    -1/3*(sin(phi)*R-(-R^2*cos(phi)^2+a^2)^(1/2))^3

P  :=  int(rho^2,
           rho = sin(phi)*R - sqrt(sin(phi)^2*R^2 - R^2 + a^2)
                 ..
                 sin(phi)*R + sqrt(sin(phi)^2*R^2 - R^2 + a^2));

     1/3*(sin(phi)*R+(sin(phi)^2*R^2-R^2+a^2)^(1/2))^3
    -1/3*(sin(phi)*R-(sin(phi)^2*R^2-R^2+a^2)^(1/2))^3

# Those are the same. The problem is with the outer integral.
simplify(PP-P);

                     0

int(P*sin(phi),phi=arccos(a/R)..Pi - arccos(a/R), method=elliptictrig);

                  1       2
                - - R Pi a 
                  2        

# The following returns unevaluated.
int(PP*sin(phi),phi=arccos(a/R)..Pi - arccos(a/R), method=elliptictrig):
op(0,%);

                    int

int(P*sin(phi),phi=arccos(a/R)..Pi - arccos(a/R), method=ftoc);

                        2
                  R Pi a 

int(PP*sin(phi),phi=arccos(a/R)..Pi - arccos(a/R), method=ftoc);

                        2
                  R Pi a 

When no method is specified the integration of ("simplified" form) P produces the wrong result from the ellitictric method, which gets selected and returned automatically. When no method is specified the integration of the (unsimplified form) PP produces the acceptable result from the ftoc method, since the elliptictrig method fails to produce anything.

Thus a workaround is to force the ftoc method for such examples that otherwise go wrong with the elliptictrig method. That may not always convenient to set up, however.

I have previously seen several similar problematic examples with that method. I'll submit a bug report.

int_oddity.mw

Here is the original nested integral, using 2D Input and with the inner integral typeset and the outer integral forcing method=ftoc. The limits of integration are as in the original, and the desired result attains.

int_ftoc.mw

My personal preference is for a cutout view, so that any relationship between the particular edges of the washers with the surfaces can be seen.

In order to keep the GUI responsive, it helps considerably to keep the frame data smaller. By using a grid=[2,49] option for the washer construction they can retain the same look as default, while taking less memory resources. The benefit is great enough that a modest number of washers can be shown together in the frames.

restart;
F:=proc(yy, ea) local hh,opts,phi,r,y; uses plots;
  opts := style=surface, color=gold, grid=[2,49];
  display(seq([
    display(plot3d([[r*cos(phi),r*sin(phi),y],
                    [r*cos(phi),r*sin(phi),y+h]],
                   r=y^2..sqrt(8*y), phi=0..ea, opts),
            plot3d([[y^2*cos(phi),y^2*sin(phi),hh],
                   [sqrt(8*y)*cos(phi),sqrt(8*y)*sin(phi),hh]],
                   hh=y..y+h, phi=0..ea, opts)),
    display(polygonplot3d([[y^2,0,y],[y^2,0,y+h],
                           [sqrt(8*y),0,y+h],[sqrt(8*y),0,y]]),
            polygonplot3d([[y^2*cos(ea),y^2*sin(ea),y],
                           [y^2*cos(ea),y^2*sin(ea),y+h],
                           [sqrt(8*y)*cos(ea),sqrt(8*y)*sin(ea),y+h],
                           [sqrt(8*y)*cos(ea),sqrt(8*y)*sin(ea),y]]),
      transparency=0.5)][], y=0..yy, h));
end proc:
numframes:=15: h:=2.0/(numframes): endangle:=4*Pi/3:
P:=plot3d([[r*cos(phi),r*sin(phi),r^2/8],[r*cos(phi),r*sin(phi),r^(1/2)]],
          r=0..4, phi=0..endangle,
          style=surface, color=["Green","Blue"], transparency=0.75):
plots:-animate(F,[y, endangle], y=0..2-h, frames=numframes,
               background=P, paraminfo=false,
               lightmodel=light2, glossiness=1.0, scaling=constrained,
               axes=normal, labels=[z,x,y], orientation=[-60,70]);

washers_ac.mw

Naturally, it can be done similarly for shells.

@mikerostructure Yes, I have Maple 2020.0, and I use right-click to export 3D plots on a regular basis.

I usually export to .png format. If your problem is with .eps export specifically then this would a good time to mention it. (For me, .eps export of your simple example does make my machine's CPU ramp up for a couple of seconds -- the export file is not small. But it doesn't crash.)

Another possible point of interest may be whether you have a particular language setting in Maple's GUI.

In your followup comments you've started mentioning a problematic file. Is there a reason that you cannot upload and attach that particular file to a Comment here, using the green up-arrow of the Mapleprimes editor? If the problem occurs when you re-execute it from scratch then perhaps you could even "Edit/Remove All Output" within it, for an easier, smaller upload.

This is numerical error, yes. The particular behavior depends on the working precision.

Your Document has a line that assigns to the name Digit . If you changes that to be Digits instead then you can easily adjust the working precision and see what happens when it is run with Digits=15,16,17,18, etc.

If this expression is assigned to name expr, say, then you can compare with expand(expr) and combine(expand(expr)). Compare with the values that occur in the numerator of your original expression.

Note that by default this plotting call will try and utilize evalhf to do the plotting in hardware double precision. That switches to so-called software floating-point when Digits>15.

This particular plotting command will also switch from evalhf  to sofware floats as a fallback if it detects a Float(undefined) result. There is a small range from about t=14880.0 to t=14902.7 in which the expression evaluates to 0.0 because roundoff error has occured but yet not underflow/overflow. For t above that range it will fall back to software float mode (and the curve again appears even at 0.45). This attachment demonstrates this behaviour 1_ac2.mw

The particular fashion in which numerical error occurs for this example is different for forced software floats at Digits<15 versus hardware double precision. There are ways to force software floats for Digits<=15 . (see attachment)

(Inlining documents to Mapleprimes seems to be broken right now.) Here is an attachment that demonstrates that forcing of software floats at lower Digits : Download 1_ac.mw

If your plot was constructed using the plot command then all you have to do is look at the product documentation, in order to see how to specify a particular range for the independent variable.

For example the Help page for the plot command shows several examples of this basic functionality of specifying the range over which to plot.

Or, was your plot generated in some other way than typing in a command? Did you use an item from the context-menu?

Just for fun, making it look closed.

I did this using Maple 16. If your version is older and this doesn't work then tell us what version you use.

F:=plottools[transform]((x,y)->[x,0,y+x^2/8])(plot(r^(1/2)-r^2/8,r=0..4,
                                                   filled,transparency=0.0)):
P := A -> plots:-display(orientation=[-45,70],
            plot3d([[r,theta,r^(1/2)],[r,theta,r^2/8]],
                   r=0..4, theta=0..A, coords=cylindrical,
                   orientation=[-50,80], lightmodel=light1,
                   grid=[40,max(2,ceil(evalf(49*A/(2*Pi))))],
                   scaling=constrained, axes=normal, style=patch),
            `if`(A<6.28,[F,plottools[rotate](F,A,[[0,0,0],[0,0,1]])][],NULL)):
plots[animate](P, [A], A=0..2*Pi, frames=17, paraminfo=false);

Download solid_rev.mw

@Mike Graber You wrote this:

        top level.

        A:=7:

        B:=A:

        A:=9:

        C:=A:

        At this point in the program, both B and C have the value 9. 

That code, by itself, will result in B having 7 as its assigned value (ie. B will evaluate to 7), and not 9.

Perhaps you've described what you want to see happen, and not what you've seen happen. If that is what you want then I think you could have expressed it more clearly.

In that case what you want is for B to get assigned a reference to A by name (and not its prior value), so that it subsequently evaluates to whatever value A subsequently evaluates (in an ongoing manner). Here is an example of that,

A := 7;

           A := 7

B := 'A';

           B := A

A := 9;

           A := 9

B;

              9

Your use of the phrase "by address" could be better replaced with "by reference".

To access or refer to a variable by its name -- rather than by its value -- use single right-quotes. Those are sometimes called unevaluation quotes in Maple context, as they delay the evaluation of the name.

See Chapter 3 of the Programming Manual (esp. Sections 3.3 and 3.4) as one place to get some background.

See also the Help pages for Topics procedure and evaluation for some additional information on specifc behavior.

These are not dumb questions. But you really ought to provide the full details of the preliminary computations.

Use the green up-arrow in the Mapleprimes editor to upload and attach your actual worksheet.

Use colon-equals := to make an assignment to a name. Your original image (now edited in the original Question) showed phi = 10 rather than an assignment statement.

By using only equals = you create an equation, which does not assign any value to the name.

Ie ,

phi:=10

versus

phi=10

The second of those -- by itself--  has not affect on later use of phi.

You can, however, utilize the equation phi=10 within an `eval` call.

As for `assume` it is tricky. Best if you show your problem case in full. My guess is that you are running into the fact that assumptions do not prevent any assignment, and that placing new assumptions (using `assume`) and assignment will wipe previous `assume` effects. See the `additionally` command. Or utilize `assuming` instead.

You already knew about assignment statements, and I suspect that you also knew about ways to concatenate names. You've put these actions together.

Note that it is only practical in a modest number, and since the generated names are always global it's far less useful as a programming technique.

Also, each global new name will get checked against the assigned names stored in all .mla and .m files in libname (once at least, up to first evaluation). That seems to scale roughly linearly, IIRC. I once saw an example with a 15sec wait in the computation where the kernel was doing nothing but i/o overhead.

So it's sorta, kinda ok as an off-the-cuff thing at the top-level, but in my opinion doesn't belong in a group of solid programming techniques for serious production line work.

You might consider this (although nothing is fool-proof),

  select(u->is(u,real), evalc~(allvalues~(coef7))) assuming alpha[1,2]::real, alpha[2,6]::real;

or, since it is only the two assumed names as dependent names,

  select(u->is(u,real), evalc~(allvalues~(coef7))) assuming real;

I declared z as local outside the do-loop as, alas, the 2D Input parser doesn't accept your original.

Download G2_ac.mw

You might consider editing your procedures in 1D plaintext Maple Notation (in Execution Groups, or Code Edit Regions).

ps. I am curious, is your task going to go anywhere near the LinearAlgebra[Generic] subpackage? (I know nothing of "Domineering".)

You haven't stated why LinearAlgebra:-LinearSolve is not adequate for your computations. You haven't mentioned whether you are doing purely floating-point computations, and if so whether those are purely real.

I would not be surprised if you needed some large examples (or a great many of them) before the LinearSolve overhead became a bigger bottleneck than anything else in your code -- as a tentative judgement based from general experience with member's code.

But, below is an example that uses the floating-point real tridiagonal linear solver from LAPACK.

This uses two functions, one to do the LU decomposition and another to do the back/forward substitutions for one or more RHSs.

Hopefully I got the integer array's widths correct -- I might have forgetten some detail of the linkage model in use.

As I mention below, this has a little extra overhead (probably avoidable, using other access means).

This was run on Linux. To run on MS-Windows you'd have to changed the name of the external shared library, from "libmkl_rt.so" to the appropriate .dll. And even then it might run amok of MS-Windows awkwardness.

Download DGTTRF.mw