[MITgcm-support] Icefront returnin NaN based on topology.

Tue Feb 18 14:24:00 EST 2014

Hi Martin

Attached is my gendata.m code which comes up with the topology, which is actually read from a CSV. The size of the simulation is small, 120 x 48 in the xy plane with grid size of 50m, and 10 intervals of 15.5 m in the z direction.

Also attached are the major data files. Other files such as the initial temperature and surface heat flux have been tested a lot and I believe are not the issue here. It is something in the icefront package which means it works OK when I turn it off, or if I turn it on and set the icelength input file to zero everywhere, but if I make this non-zero, as is in the gendata.m attached, it breaks down.

I've also not got any initial results to pull a pressure distribution from, so as my lake is very slowly moving, I have just used pressure = rho*g*h.
In Size.h if I replace Nr everywhere with Nz then the model returns no output. (Or at least, none of the T.000000050.001.001.data files for T, U V W ....) but from looking at tutorial_Deep_convection I thought it would be ok to leave this as Nr seeing as it is done in that tutorial, which uses a similarly sized grid and Cartesian coordinates. 

Thank you so much for any assistance!

Harry

Here is the contents of the gendata.m file attached:

close all
ieee='b';
accuracy='real*8';

nx=120;
ny=48;
nz=10;
dz=15.5;
dy=50;
dx=dy;
x=csvread('120by48topology.csv'); %This topology is between 0 and 150 meters, with zero's all around the outside so all boundaries are fixed.

%So making our topology...
h=-1*x;
h=2*round(h/2); %Making h multiples of 2 rather than rational numbers.
h=h';
%And the overall topology
fid=fopen('topog.box','w',ieee); fwrite(fid,h,accuracy); fclose(fid);

%Now, making our shelfice input file
icet= abs(h)>0;
icetopo=icet*0.0; 
icetopo(1:55,1:48)=1*h(1:55,1:48); %So now we make the shelf ice the same depth as the topology over the top 55 spaces in domain
fid=fopen('icetopo.exp1','w','b'); fwrite(fid,icetopo,accuracy);fclose(fid);

%so, now we generate our icefront.We just set it to be equal to the topology height at the front of the shelfice file, as is done in isomip
iceheight=0*h;
iceheight(56,1:48)=-1*h(55,1:48);%%%%%%%%% 
fid=fopen('ifront.xuyun','w','b'); fwrite(fid,iceheight,accuracy); fclose(fid);

%Now making the icelength file
icelength=iceheight*0.0+(dz/dy)*(iceheight~=0); %So, it will only have data in the same places that the iceheight does
fid=fopen('ilength.xuyun','w','b'); fwrite(fid,icelength,accuracy); fclose(fid);

%Now making the potential distribution using simple pressure = rho*9.81*height
potential=ones(nx,ny,nz);
for k=1:nz
potential(:,:,k)=9.81*dz*998*potential(:,:,k)*(k-0.5);
end
fid=fopen('potential.jmd95z','w','b'); fwrite(fid,potential,accuracy);fclose(fid);

> From: Martin.Losch at awi.de
> Date: Tue, 18 Feb 2014 13:17:27 +0100
> To: mitgcm-support at mitgcm.org
> Subject: Re: [MITgcm-support] Icefront returnin NaN based on topology.
> 
> Harry,
> 
> I think it would be useful to know about the details of the configuration. It would also help to have the illustrative example that you describe as a “code” and “input” directory (I am assuming it’s a very small domain), including the matlab (?) file to generate the different input files.
> 
> Martin
> 
> On Feb 16, 2014, at 2:41 AM, Harry A <gmtonys3 at hotmail.com> wrote:
> 
> > Hi Everyone
> > 
> > I'm working on a model including a glacier, so I've been implementing the Icefront package to add this functionality.
> > However, when I am setting up the topology of my submerged ice, I've had NaN results for all outputted velocities and temperatures turning up at strange times.
> > 
> > I have a bottom topology over my whole domain with closed boundaries (Depth = 0 all around) To begin with, I just made a 3x3 grid in the middle of my domain ice, but depending on where I place it the model with either return all NaN or return good results. It seems even if I make this ice shallow (say 10 meters deep in 100m water) or as deep as the water at the certain location, or even deeper, I still get failure with some locations and not others.
> > 
> > To highlight the problem, one example is I can make a single grid point in the middle of my domain have ice 70.4952 meters deep (with the bathymetry file being exactly 70m deep in that location) and it will work, however if I make this 70.4953 or 70.4951 the model returns all NaNs.
> > 
> > Does anyone have any idea what could be causing these issues? Am I correct in saying the Ice depth file for Icefront ( ICEFRONTdepthFile) meant to be positive values ranging from 0 up to the depth of water at any given grid point ? And the ICEFRONTlengthFile should be the facial area of a grid element in the horizontal direction, so If I'm using dX=50, dY=50, dZ = 15, it should be 15*50 = 750?
> > 
> > Thanks for any assistance, I've been running through matlab scrips generating ice topologies for 6 hours now and there seems to be almost no rhyme or reason as to why some situations fail and some work just fine...
> > 
> > Harry
> > _______________________________________________
> > MITgcm-support mailing list
> > MITgcm-support at mitgcm.org
> > http://mitgcm.org/mailman/listinfo/mitgcm-support
> 
> 
> _______________________________________________
> MITgcm-support mailing list
> MITgcm-support at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-support

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140218/a3cc580c/attachment-0001.htm>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: data
Type: application/octet-stream
Size: 2776 bytes
Desc: not available
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140218/a3cc580c/attachment-0003.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: data.icefront
Type: application/octet-stream
Size: 333 bytes
Desc: not available
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140218/a3cc580c/attachment-0004.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: data.pkg
Type: application/x-newton-compatible-pkg
Size: 71 bytes
Desc: not available
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140218/a3cc580c/attachment-0001.bin>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: data.shelfice
Type: application/octet-stream
Size: 643 bytes
Desc: not available
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140218/a3cc580c/attachment-0005.obj>
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: CPP_OPTIONS.h
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140218/a3cc580c/attachment-0002.h>
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: packages.conf
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140218/a3cc580c/attachment-0002.el>
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: SIZE.h
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140218/a3cc580c/attachment-0003.h>
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: gendata.m
URL: <http://mitgcm.org/pipermail/mitgcm-support/attachments/20140218/a3cc580c/attachment-0003.el>