Merge pull request #8 from gaelforget/AbstractArrayAsSupertype

Abstract array as supertype

Author: gaelforget

Project.toml

@@ -1,7 +1,7 @@
 name = "MeshArrays"
 uuid = "28126448-9d44-11e8-08f5-e59c47736fcb"
 authors = ["gaelforget <[email protected]>"]
-version = "0.1.4"
+version = "0.1.5"
 
 [deps]
 Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"

src/MeshArrays.jl

@@ -12,7 +12,8 @@ include("gcmfaces_convert.jl");
 include("gcmfaces_IO.jl");
 include("gcmfaces_demo.jl");
 
-export gcmfaces, exchange, gradient, convergence, smooth, mask, fsize
+export AbstractGcmfaces, gcmfaces, gcmsubset, fsize, fijind
+export exchange, gradient, convergence, smooth, mask
 export GCMGridSpec, GCMGridLoad, GCMGridOnes
 export demo1, demo2
 #The following functions rely on grid specs; currently via global vars.

src/gcmfaces_calc.jl

@@ -9,8 +9,8 @@ end
 function gradient(inFLD::gcmfaces,doDIV::Bool)
 
 exFLD=exchange(inFLD,1)
-dFLDdx=gcmfaces(inFLD.nFaces,inFLD.grTopo)
-dFLDdy=gcmfaces(inFLD.nFaces,inFLD.grTopo)
+dFLDdx=similar(inFLD)
+dFLDdy=similar(inFLD)
 
 for a=1:inFLD.nFaces;
   (s1,s2)=fsize(exFLD,a)
@@ -32,8 +32,8 @@ end
 function gradient(inFLD::gcmfaces,iDXC::gcmfaces,iDYC::gcmfaces)
 
 exFLD=exchange(inFLD,1)
-dFLDdx=gcmfaces(inFLD.nFaces,inFLD.grTopo)
-dFLDdy=gcmfaces(inFLD.nFaces,inFLD.grTopo)
+dFLDdx=similar(inFLD)
+dFLDdy=similar(inFLD)
 
 for a=1:inFLD.nFaces;
   (s1,s2)=fsize(exFLD,a)
@@ -55,7 +55,7 @@ return fldmsk
 end
 
 function mask(fld::gcmfaces, val::Number)
-  fldmsk=gcmfaces(fld.nFaces,fld.grTopo)
+  fldmsk=similar(fld)
   for a=1:fld.nFaces
     tmp1=copy(fld.f[a])
     replace!(x -> !isfinite(x) ? val : x, tmp1 )
@@ -65,7 +65,7 @@ function mask(fld::gcmfaces, val::Number)
 end
 
 function mask(fld::gcmfaces, val::Number, noval::Number)
-  fldmsk=gcmfaces(fld.nFaces,fld.grTopo)
+  fldmsk=similar(fld)
   for a=1:fld.nFaces
     tmp1=copy(fld.f[a])
     replace!(x -> x==noval ? val : x, tmp1  )
@@ -83,7 +83,7 @@ function convergence(uFLD::gcmfaces,vFLD::gcmfaces);
 #  if otherwise then something this may be needed:
 #  uFLD=mask(uFLD,0.0); vFLD=mask(vFLD,0.0);
 
-CONV=gcmfaces(uFLD.nFaces,uFLD.grTopo)
+CONV=similar(uFLD)
 
 (tmpU,tmpV)=exch_UV(uFLD,vFLD)
 for a=1:tmpU.nFaces
@@ -118,8 +118,8 @@ mskW=mask(mskW,0.0)
 mskS=mask(mskS,0.0)
 
 #get inverse grid spacing:
-iDXC=gcmfaces(FLD.nFaces,FLD.grTopo)
-iDYC=gcmfaces(FLD.nFaces,FLD.grTopo)
+iDXC=similar(FLD)
+iDYC=similar(FLD)
 for a=1:FLD.nFaces;
   iDXC.f[a]=1.0./MeshArrays.DXC.f[a]
   iDYC.f[a]=1.0./MeshArrays.DYC.f[a]
@@ -147,8 +147,8 @@ dtFac=dt*mskC/MeshArrays.RAC;
 #loop:
 for it=1:nbt
   (dTdxAtU,dTdyAtV)=gradient(FLD,iDXC,iDYC);
-  tmpU=gcmfaces(FLD.nFaces,FLD.grTopo)
-  tmpV=gcmfaces(FLD.nFaces,FLD.grTopo)
+  tmpU=similar(FLD)
+  tmpV=similar(FLD)
   for a=1:FLD.nFaces
       tmpU.f[a]=dTdxAtU.f[a].*KuxFac.f[a];
       tmpV.f[a]=dTdyAtV.f[a].*KvyFac.f[a];

src/gcmfaces_convert.jl

@@ -23,22 +23,23 @@ return arr;
 end
 
 ## convert2array ##
-function convert2array(fld::Array)
+function convert2array(fld::Array{T,N}) where {T,N}
 
 grTopo=MeshArrays.grTopo;
 nFaces=MeshArrays.nFaces;
 
-v1=Array{Any,1}(undef,nFaces);
+v1=Array{Array{T,N},1}(undef,nFaces);
+N>2 ? error("N>2 case not implemented yet") : nothing
 
 if MeshArrays.grTopo=="llc";
-    (n1,n2)=MeshArrays.facesSize[1,:];
+    (n1,n2)=MeshArrays.facesSize[1];
     v1[1]=fld[1:n1,1:n2];
     v1[2]=fld[n1+1:n1*2,1:n2];
     v1[3]=rotr90(fld[1:n1,n2+1:n2+n1]);
     v1[4]=rotl90(fld[n1*2+1:n1*3,1:n2]);
     v1[5]=rotl90(fld[n1*3+1:n1*4,1:n2]);
 elseif MeshArrays.grTopo=="cs";
-    (n1,n2)=MeshArrays.facesSize[1,:];
+    (n1,n2)=MeshArrays.facesSize[1];
     v1[1]=fld[1:n1,n1+1:n1+n2];
     v1[2]=fld[n1+1:2*n1,n1+1:n1+n2];
     v1[3]=rotr90(fld[1:n1,n1+n2+1:n2+n1*2]);
@@ -66,7 +67,7 @@ function convert2gcmfaces(fld::gcmfaces)
 
     aa=0;bb=0;
     for iFace=1:nFaces;
-        aa=aa+prod(facesSize[iFace,:]);
+        aa=aa+prod(facesSize[iFace]);
         bb=bb+prod(size(fld.f[iFace]));
     end;
     n3=Int64(bb/aa);
@@ -75,8 +76,8 @@ function convert2gcmfaces(fld::gcmfaces)
     i0=0; i1=0;
     for iFace=1:nFaces;
         i0=i1+1;
-        nn=facesSize[iFace,1];
-        mm=facesSize[iFace,2];
+        nn=facesSize[iFace][1];
+        mm=facesSize[iFace][2];
         i1=i1+nn*mm;
         v11[i0:i1,:]=reshape(fld.f[iFace],(nn*mm,n3));
     end;
@@ -103,7 +104,7 @@ v1=Array{Array{eltype(fld),ndims(fld)}}(undef,nFaces);
 i0=0; i1=0;
 for iFace=1:nFaces
   i0=i1+1;
-  nn=facesSize[iFace,1]; mm=facesSize[iFace,2];
+  nn=facesSize[iFace][1]; mm=facesSize[iFace][2];
   i1=i1+nn*mm;
   if n3>1;
     v1[iFace]=reshape(v0[i0:i1,:],(nn,mm,n3));

src/gcmfaces_exch.jl

@@ -81,7 +81,7 @@ fillval=0.0
 #step 1
 
 s=size.(fld.f);
-FLD=gcmfaces(1,"ll");
+FLD=similar(fld);
 FLD.f[1]=fill(fillval,s[1].+2N);
 @views FLD.f[1][N+1:N+s[1][1],N+1:N+s[1][2]]=fld.f[1];
 
@@ -118,8 +118,8 @@ fillval=0.0
 #step 1
 
 s=size.(fldU.f);
-FLDU=gcmfaces(1,"ll");
-FLDV=gcmfaces(1,"ll");
+FLDU=similar(fldU);
+FLDV=similar(fldV);
 
 FLDU.f[1]=fill(fillval,s[1][1]+1,s[1][2]);
 FLDV.f[1]=fill(fillval,s[1][1],s[1][2]+1);
@@ -148,7 +148,7 @@ s=size.(fld.f);
 nf=fld.nFaces;
 nf==5 ? s=vcat(s,s[3]) : nothing
 tp=fld.grTopo;
-FLD=gcmfaces(nf,tp);
+FLD=similar(fld);
 
 for i=1:nf; FLD.f[i]=fill(fillval,s[i].+2N); end;
 #code below yields strange, seemingly incorrect results:
@@ -196,8 +196,8 @@ s=size.(fldU.f);
 nf=fldU.nFaces;
 nf==5 ? s=vcat(s,s[3]) : nothing
 tp=fldU.grTopo;
-FLDU=gcmfaces(nf,tp);
-FLDV=gcmfaces(nf,tp);
+FLDU=similar(fldU);
+FLDV=similar(fldV);
 
 for i=1:nf;
  FLDU.f[i]=fill(fillval,s[i].+2N);
@@ -250,8 +250,8 @@ s=size.(fldU.f);
 nf=fldU.nFaces;
 nf==5 ? s=vcat(s,s[3]) : nothing
 tp=fldU.grTopo;
-FLDU=gcmfaces(nf,tp);
-FLDV=gcmfaces(nf,tp);
+FLDU=similar(fldU);
+FLDV=similar(fldV);
 
 for i=1:nf
   FLDU.f[i]=fill(fillval,s[i][1]+1,s[i][2]);

src/gcmfaces_grids.jl

@@ -14,21 +14,21 @@ if gridName=="LLC90";
     nFaces=5;
     grTopo="llc";
     ioSize=[90 1170];
-    facesSize=[[90 270];[90 270];[90 90];[270 90];[270 90]];
+    facesSize=[(90, 270), (90, 270), (90, 90), (270, 90), (270, 90)]
     ioPrec=Float64;
 elseif gridName=="CS32";
     grDir="GRID_CS32/";
     nFaces=6;
     grTopo="cs";
     ioSize=[32 192];
-    facesSize=[[32 32];[32 32];[32 32];[32 32];[32 32];[32 32]];
+    facesSize=[(32, 32), (32, 32), (32, 32), (32, 32), (32, 32), (32, 32)]
     ioPrec=Float32;
 elseif gridName=="LL360";
     grDir="GRID_LL360/";
     nFaces=1;
     grTopo="ll";
     ioSize=[360 160];
-    facesSize=[360 160];
+    facesSize=[(360, 160)]
     ioPrec=Float32;
 else;
     error("unknown gridName case");
@@ -63,9 +63,10 @@ function GCMGridOnes(grTp,nF,nP)
     grDir=""
     grTopo=grTp
     nFaces=nF
-    ioSize=[nP nP*nF];
-    facesSize=fill(nP,(nF,2));
-    ioPrec=Float32;
+    ioSize=[nP nP*nF]
+    facesSize=Array{NTuple{2, Int},1}(undef,nFaces)
+    facesSize[:].=[(nP,nP)]
+    ioPrec=Float32
 
     list0=("XC","XG","YC","YG","RAC","RAZ","DXC","DXG","DYC","DYG","hFacC","hFacS","hFacW","Depth");
     for ii=1:length(list0);