Question: How to create 3Dplot of 4dim. function using 3Dplot/structure

Hello,

I have problem with creating 3D plot of this function:

a:=(alpha,H,gamma)->(Pi^2/(4*alpha)-Pi/2+1)*2/H*(sin(gamma/2)^2):

where

alpha, H, gamma > 0 ... independent variables,

a ... dependent variable.

3D plot should be an set of 3D points in "alpha", "H" and "gamma" range (x, y, z coord.), value of dependent variable "a" should be indicated by color of each point.

I tried to use command PLOT3D - ISOSURFACE, but it still doesnt work because of some error in plot command and I dont know what is wrong... -

please see my code:

 #####################################

a:=(alpha,H,gamma)->(Pi^2/(4*alpha)-Pi/2+1)*2/H*(sin(gamma/2)^2):

D2R:=evalf(Pi/180);

# alpha [rad] # range and step of indep. var. "alpha"
 imin:=16*D2R;
 delta_i:=2*D2R;
 imax:=30*D2R;
 
 # H [m] # range and step of indep. var. "H"
 jmin:=0.08;
 delta_j:=0.005;
 jmax:=0.15;
 
 # gamma [rad] # range and step of indep. var. "gamma"
 kmin:=30*D2R;
 delta_k:=2*D2R;
 kmax:=40*D2R;

# max. values of variables in repetition statement 'for' #

mmax:=trunc((imax-imin)/delta_i+1);
nmax:=trunc((jmax-jmin)/delta_j+1);
pmax:=trunc((kmax-kmin)/delta_k+1);
 
# repetition statement for definition of all 4 values of each 3D point: #
 for m from 1 to mmax do
   for n from 1 to nmax do
     for p from 1 to pmax do

       i:=imin+(m-1)*delta_i;
       j:=jmin+(n-1)*delta_j;
       k:=kmin+(p-1)*delta_k;
     
       A[m,n,p]:=evalf([i,j,k,a(i,j,k)]);
         
     end do;
   end do;
 end do;

 
 unassign('i','j','k','m','n','p');


# This should be an 'list of m lists. Each sublist in turn contains n lists with p elements, each of which is a list [xijk, yijk, zijk, fijk],

representing the (x, y, z) coordinates and the function value of grid point (i, j, k)' - for more see help for plot3d/structure - ISOSURFACE: #

array_mnp:=[[[A[m,n,p] $ p=1..pmax] $ n=1..nmax] $ m=1..mmax];

PLOT3D(ISOSURFACE(array_mnp);

 #####################################

Using this plot command generate only "Plotting error, empty plot" though all the values of "array_mnp" are float-type.

------------------------------------------------------------------------

My second version of the code use an array AA (against list "array_mnp" in the 1st version) in plot command:

 #####################################

AA:=hfarray(1..mmax,1..nmax,1..pmax,1..mmax*nmax*pmax); # Here is probably the problem, I'm not sure with this declaration #

for m from 1 to mmax do
  for n from 1 to nmax do
     for p from 1 to pmax do
 
       i:=imin+(m-1)*delta_i;
       j:=jmin+(n-1)*delta_j;
       k:=kmin+(p-1)*delta_k;
 
       AA[m,n,p,1]:=i;
       AA[m,n,p,2]:=j;
       AA[m,n,p,3]:=k;
 
       AA[m,n,p,4]:=evalf(a(i,j,k));
 
     end do;
   end do;
 end do;

unassign('i','j','k','m','n','p');

PLOT3D(ISOSURFACE(AA));

 #####################################

The only result of this second code is "Plotting error, ISOSURFACE Arrays must be NxMxLx4". I dont know what is wrong again...

Please help me with my problem:

1) what is wrong in the first version code,

or

2) if I should only correct declaration of array "AA"

or

3) if I have to use other command then PLOT3D - ISOSURFACE


.... Thank you a lot!!


O.M.