ucftools/matlab/write_uvwp_newmesh_ucf.m

134 lines
5.6 KiB
Matlab

function [] = write_uvwp_newmesh_ucf(file,...
u,v,w,p,...
xuo,yuo,zuo,xvo,yvo,zvo,xwo,ywo,zwo,xpo,ypo,zpo,...
xun,yun,zun,xvn,yvn,zvn,xwn,ywn,zwn,xpn,ypn,zpn,...
ibegu,jbegu,kbegu,iendu,jendu,kendu,...
ibegv,jbegv,kbegv,iendv,jendv,kendv,...
ibegw,jbegw,kbegw,iendw,jendw,kendw,...
ibegp,jbegp,kbegp,iendp,jendp,kendp,...
time,varargin)
% [] = write_uvwp_newmesh_ucf(file,...
% u,v,w,p,...
% xuo,yuo,zuo,xvo,yvo,zvo,xwo,ywo,zwo,xpo,ypo,zpo,...
% xun,yun,zun,xvn,yvn,zvn,xwn,ywn,zwn,xpn,ypn,zpn,...
% ibegu,jbegu,kbegu,iendu,jendu,kendu,...
% ibegv,jbegv,kbegv,iendv,jendv,kendv,...
% ibegw,jbegw,kbegw,iendw,jendw,kendw,...
% ibegp,jbegp,kbegp,iendp,jendp,kendp,...
% time,varargin)
% Writes a complete uvwp field to UCF chunks with a different mesh.
% Caution: extrapolation might be performed (periodicity is not checked)
% Input
% file output file
% u,v,w,p complete uvwp fields on old grid
% xuo,yuo,zuo,... old mesh
% xun,yun,zun,... new mesh
% ibegu,jbegu,kbegu,... new processor grid
% time simulation time
% Optional input (key/value pair)
% method interpolation method (default: 'cubic')
% endian endianess of the file as used by fopen (default: 'n')
% verbosity verbose output? (default: no)
% debug debug output? (default: no)
% Parse variable input
par = inputParser;
addParamValue(par,'method','cubic',@ischar);
addParamValue(par,'endian','n',@ischar);
addParamValue(par,'verbosity',0,@isnumeric);
addParamValue(par,'debug',0,@isnumeric);
parse(par,varargin{:});
% Parse filename
[fdir,fbase,fext] = fileparts(file);
fbasepath = fullfile(fdir,fbase);
nxprocs = numel(ibegu);
nyprocs = numel(jbegu);
nzprocs = numel(kbegu);
nxu = numel(xun); nxv = numel(xvn); nxw = numel(xwn); nxp = numel(xpn);
nyu = numel(yun); nyv = numel(yvn); nyw = numel(ywn); nyp = numel(ypn);
nzu = numel(zun); nzv = numel(zvn); nzw = numel(zwn); nzp = numel(zpn);
% Create interpolants
giu = griddedInterpolant({xuo,yuo,zuo},u,par.Results.method);
giv = griddedInterpolant({xvo,yvo,zvo},v,par.Results.method);
giw = griddedInterpolant({xwo,ywo,zwo},w,par.Results.method);
gip = griddedInterpolant({xpo,ypo,zpo},p,par.Results.method);
for ixproc=0:nxprocs-1
for iyproc=0:nyprocs-1
for izproc=0:nzprocs-1
[iproc] = locfun_proc_id(ixproc,iyproc,izproc,nxprocs,nyprocs,nzprocs);
fname = sprintf('%s.%05d',fbasepath,iproc);
% Get processor boundaries
ibu = ibegu(ixproc+1); ibv = ibegv(ixproc+1); ibw = ibegw(ixproc+1); ibp = ibegp(ixproc+1);
ieu = iendu(ixproc+1); iev = iendv(ixproc+1); iew = iendw(ixproc+1); iep = iendp(ixproc+1);
jbu = jbegu(iyproc+1); jbv = jbegv(iyproc+1); jbw = jbegw(iyproc+1); jbp = jbegp(iyproc+1);
jeu = jendu(iyproc+1); jev = jendv(iyproc+1); jew = jendw(iyproc+1); jep = jendp(iyproc+1);
kbu = kbegu(izproc+1); kbv = kbegv(izproc+1); kbw = kbegw(izproc+1); kbp = kbegp(izproc+1);
keu = kendu(izproc+1); kev = kendv(izproc+1); kew = kendw(izproc+1); kep = kendp(izproc+1);
% Create output file
obj = ucf('verbosity',par.Results.verbosity,'debug',par.Results.debug);
obj.create(fname,'type','field','rank',[iproc,ixproc,iyproc,izproc],'endian',par.Results.endian);
obj.appendStep(time);
ighost = 0;
% Interpolate to refined chunk and write
data = giu({xun(ibu:ieu),yun(jbu:jeu),zun(kbu:keu)});
if any(isnan(data))
warning('NaN detected: u, proc #%d, %d occurences',iproc,sum(isnan(data(:))));
end
nxul = size(data,1)-2*ighost;
nyul = size(data,2)-2*ighost;
nzul = size(data,3)-2*ighost;
params = int64([ighost,ibu,jbu,kbu,nxul,nyul,nzul,nxu,nyu,nzu]);
obj.appendSet(data,params);
% v
data = giv({xvn(ibv:iev),yvn(jbv:jev),zvn(kbv:kev)});
if any(isnan(data))
warning('NaN detected: v, proc #%d, %d occurences',iproc,sum(isnan(data(:))));
end
nxvl = size(data,1)-2*ighost;
nyvl = size(data,2)-2*ighost;
nzvl = size(data,3)-2*ighost;
params = int64([ighost,ibv,jbv,kbv,nxvl,nyvl,nzvl,nxv,nyv,nzv]);
obj.appendSet(data,params);
% w
data = giw({xwn(ibw:iew),ywn(jbw:jew),zwn(kbw:kew)});
if any(isnan(data))
warning('NaN detected: w, proc #%d, %d occurences',iproc,sum(isnan(data(:))));
end
nxwl = size(data,1)-2*ighost;
nywl = size(data,2)-2*ighost;
nzwl = size(data,3)-2*ighost;
params = int64([ighost,ibw,jbw,kbw,nxwl,nywl,nzwl,nxw,nyw,nzw]);
obj.appendSet(data,params);
% p
data = gip({xpn(ibp:iep),ypn(jbp:jep),zpn(kbp:kep)});
if any(isnan(data))
warning('NaN detected: p, proc #%d, %d occurences',iproc,sum(isnan(data(:))));
end
nxpl = size(data,1)-2*ighost;
nypl = size(data,2)-2*ighost;
nzpl = size(data,3)-2*ighost;
params = int64([ighost,ibp,jbp,kbp,nxpl,nypl,nzpl,nxp,nyp,nzp]);
obj.appendSet(data,params);
obj.close();
end
end
end
end
function [ind] = locfun_proc_id(ii,jj,kk,nx,ny,nz)
% local version of 'sub2ind_zero_row'
ind=ii*ny*nz+jj*nz+kk;
end