@planetmknzm 

1. Please explain in words what the code does.

2. Is it just a tangent plane and TNB frame moving along a space curve embedded on a sphere?

3. If you want me to help, first make a substantially simpler version of the same thing, 1/4 the current code size or less.

4. I have no interest in improving some lengthy code AND THEN trying to make those improvements work in Maple 12 when the earliest Maple that I have access to is Maple 2019. You'll need to get a more-recent version for me to work on it.

5. Don't email me; use MaplePrimes.

@Jno Here's another version, nearly twice as fast (not that the previous version was slow by any means), and more simply coded. I originally used ArrayTools:-AddAlongDimension because I mistakenly thought that all the commands is ArrayTools used external compiled code for efficiency. But AddAlongDimension is a 57-line Maple procedure, which is quite surprising because I replaced it with a 1-line procedure using only top-level built-in commands:

AAD:= (A,d)-> local k; add(A[(..)$d-1, k], k= [rtable_dims](A)[d])

It seemed worthwhile to make the change because that one-liner is so glaringly simple. So, here's the new code:

ANOVA_3D:= (X::And(rtable, 3 &under rtable_num_dims))-> 
local
    N:= numelems(X),
    #replacement for ArrayTools:-AddAlongDimension:
    AAD:= (A,d)-> local k; add(A[(..)$d-1, k], k= [rtable_dims](A)[d]),
    SS:= A-> add(A^~2)*numelems(A)/N,
    Xij:= AAD(X, 3),  Xi:= AAD(Xij, 2),  
    all:= SS(<add(Xi)>), SSi:= SS(Xi),  SSj:= SS(AAD(Xij, 1))    
; 
    Record(
        ("SS__P", "SS__A", "TSS")=~ (SSi, SSj, SS(X)) -~ all,
        "SS__AP"= SS(Xij) - SSi - SSj + all
    )
:

The example usage, calling protocol, etc., are identical.

From anywhere in a Maple worksheet, type Ctrl-J or Ctrl-K. This should give you a prompt for 1D input:

> 

If you copy-and-paste my code to the space in front of that prompt and press Enter, it should work.

Regarding those vectors of measurements that you mentioned earlier: Do they all have the same number of measurements? If they don't, then I think that you will need to use a (2-D) matrix of vectors like you were mentioning. It's not a big deal to make that change. 
 

The summation limits Sum(..., n= 1...n) should also be corrected.

@Teep No problem. You can certainly cut and paste the commands to a new worksheet and be done with the display problems. Those problems have nothing to do with the code. You may be able to also cut and paste your expository text.

@planetmknzm Okay, here's a lesson in plots and deconstruction: Almost all of Maple's plotting commands produce human-readable data structures, which are sent on to other software (called "renderers") for the actual plotting. This was surprising to me when I first learned it; I don't think that many software packages allow users to poke around in the internal structures.

The most-basic of all of Maple's deconstruction commands is op (which stands for operand or operands). The indets that I used in my Answer above is also a fundamental deconstructor, but you can learn a lot and do a lot just with op, so let's just discuss that here. 

Take the result of any plotting command, and assign it to a variable:

c:= plots:-spacecurve([(cos,sin)(t), t], t= -Pi..Pi):

Now start poking around in the data structure with op:

op(c);

You'll probably see the word CURVES, some big parentheses, and an array of numbers.

op(0, c);

You'll definitely see the word PLOT3D, which is the so-called zeroeth operand of any 3D plot data structure.

op(1, c);

In this case, you'll see the same as op(c) because this plot structure only has 1 operand (not counting the 0th). 

Digging deeper,

op([1,1], c);

That's the 1st operand of the 1st operand. In this case, that's the array of numbers.

op([1,0], c);

You get the word CURVES, which is the 0th operand of the first operand.

@Carl Love It should be noted that my technique above is not limited to selecting the first entries. Any one entry from each class can be selected. The following selects a random equation from each class:

R:= 1+~(rand@nops)~((AA2:= ([entries]@op~@ListTools:-Classify)(lhs, AA))):
AnyCombo:= ()-> op~(R(), AA2):
AnyCombo();

AnyCombo();

For 2D-Input the first line of code may need to be changed to

AA2:= ([entries]@op~@ListTools:-Classify)(lhs, AA):
R:= 1+~(rand@nops)~(AA2):
 

@vv Thank you. In the meantime, both the worksheet's original creator and I have noticed several other anomalies with the output display from this worksheet.

@Teep Indeed, there is something very odd that is specific to the output display in your worksheet. I don't know what it is; perhaps you made some setting, or perhaps it's some type of corruption. I had already posted a separate Question about it: "Bad output display from LPSolve".

@tomleslie I disagree with what you said in two ways.

First, I think that there are often good reasons for extracting the numeric data from a plot rather than generating it by other means. I do it very often. My most-common reason is to get the min and max of the plotted points in all dimensions in order to programmatically set a view that is slightly larger than those ranges. My second most-common reason is to track down the specific parameter values that have led to some plot anomaly.

Second, I think that the OP has provided enough detail in this case that the Question can be adequately Answered.

@TechnicalSupport For a table of unknown expressions, it's much more likely that Grid will work and Threads won't. Of course, if both could work,Threads is likely to be more efficient; but that's a very unlikely "if". Grid is also much easier for a newbie to understand, whereas Threads nearly requires a CompSci degree or equivalent. (This is not meant as a criticism of Threads, which is necessarily difficult to understand; rather than being difficult due to a poor implementation or documentation.)

You have a bug in your second constraint set. You refer to a vector Supply which hasn't been defined. I assume that you mean the vector RegionA. This is what I used in my Answer below.

So x is a 3D array. You've probably figured out by now that x can be declared in Maple via

x:= Array(1..n, 1..k, 1..r);

I understand squaring of individual elements of x. Indeed, the 1st part of your TSS formula is the sum of the squares of all elements of x, which can be done in Maple simply as add(x^~2) (note the tilde ~). But I don't know what you mean by the "squares" of "slices" of x, such as in your SS_P formula x[i, .., ..]^2/(k*r). I think that what you likely mean is the squares of the means (averages) of various rows, columns, planes of x. In that case, we'd need to first add the subsection of x, divide by the number of elements, then square the result.

@Tamour_Zubair I think that something can be done, but I think that you might have made a mistake in entering the system. One of Z_0(t), Z_1(t), Z_2(t) seems to be "extra". Did you really intend for the system to contain all three?

@acer Since acer couldn't upload his worksheet (due to limitations of MaplePrimes), here is a transcription:

restart;
with(Student:-Calculus1): with(DocumentTools): with(DocumentTools:-Layout):

P:=VolumeOfRevolution(2*x^2+1,x=0..1,volumeoptions=[color= ["Niagara 1"],
                      transparency=0.5],orientation=[0,90,0],labels=[z,x,y],
                      view=[-1..1,-2..2,0..3],output=plot,axis=vertical,
                      scaling=constrained):
P;

InsertContent(Worksheet(Group(Output(Textfield(InlinePlot(
   subsindets(plots:-display(P,size=[600,500]),'specfunc'(_TYPESET),
              u->_TYPESET(Typesetting:-Typeset~([op(u)])[])),
   ':-scale'=1.9,ytrans=2)))))):

@achreftabet The integral in your photo is a definite integral with respect to theta. Thus, its value cannot depend on theta.