Changeset 7810

Timestamp:

04/11/11 11:55:54 (14 years ago)

Author:

Mathieu Morlighem

Message:

Added TEST for permissions

Location:

issm/trunk

Files:

• TEST (added)
• examples/Bumps/Bump1_surface_bed/Square.par (modified) (5 diffs)
• examples/Bumps/Bump1_surface_bed/runme.m (modified) (1 diff)
• externalpackages/export_fig/ghostscript.m (modified) (1 diff)
• src/m/utils/Mesh/FixMesh.m (modified) (2 diffs)
• src/mex/Test/Test.cpp (modified) (2 diffs)

issm/trunk/examples/Bumps/Bump1_surface_bed/Square.par

@@ -13 +13 @@
 
 %Solution parameters
-        %parallelization 
-        md.cluster='none';
-        md.np=2;
-        md.time=1;
-        md.exclusive=0;
+        md.cluster=eval(oshostname());
 
         %statics
         md.lowmem=1;
-        md.eps_rel=0.01;
+        md.eps_res=0.05;
+        md.eps_rel=0.1;
         md.eps_abs=10; %m/yr
         md.penalty_offset=4;
@@ -36 +33 @@
         md.artificial_diffusivity=1;
 
-        %control
-        md.control_type={'drag'}; %'drag', 'B'
-        md.nsteps=5;
-        md.tolx=10^-4;
-        md.maxiter=20;
-        md.optscal=10;
-        md.fit='logarithmic'; %'absolute','relative','logarithmic'
-        md.meanvel=1000/md.yts; %1000 meters/year
-        md.epsvel=eps;
-
-
         disp('      creating thickness');
         hmin=500;
@@ -58 +44 @@
         %Add bumps to the ice sheet :
         pos=find(sqrt((md.x-500000).^2+(md.y-500000).^2)<100000); %for a cicular bump
-        md.bed(pos)=md.bed(pos)+150*cos(pi/200000*(sqrt((md.y(pos)-500000).^2+(md.x(pos)-500000).^2))); %for a circular bump
+        md.bed(pos)=md.bed(pos)+30*cos(pi/200000*(sqrt((md.y(pos)-500000).^2+(md.x(pos)-500000).^2))); %for a circular bump
 
         md.surface=md.bed+md.thickness; %the surface follows the bed so the bump is transmitted to the surface
@@ -68 +54 @@
         
         disp('      creating drag');
-        md.drag_type=2; %0 none 1 plastic 2 viscous
-        md.drag=200*ones(md.numberofgrids,1); %q=1.
+        md.drag_coefficient=200*ones(md.numberofgrids,1); %q=1.
         %Take care of iceshelves: no basal drag
         pos=find(md.elementoniceshelf);
-        md.drag(md.elements(pos,:))=0;
-        md.p=ones(md.numberofelements,1);
-        md.q=ones(md.numberofelements,1);
+        md.drag_coefficient(md.elements(pos,:))=0;
+        md.drag_p=ones(md.numberofelements,1);
+        md.drag_q=ones(md.numberofelements,1);
 
         disp('      creating temperature');
@@ -80 +65 @@
 
         disp('      creating flow law paramter');
-        md.B=paterson(md.observed_temperature);
-        md.n=3*ones(md.numberofelements,1);
+        md.rheology_B=paterson(md.observed_temperature);
+        md.rheology_n=3*ones(md.numberofelements,1);
 
         disp('      creating accumulation rates');
-        md.accumulation=ones(md.numberofgrids,1)/md.yts; %1m/a
-        md.melting=0*ones(md.numberofgrids,1)/md.yts; %1m/a
+        md.accumulation_rate=ones(md.numberofgrids,1)/md.yts; %1m/a
+        md.melting_rate=0*ones(md.numberofgrids,1)/md.yts; %1m/a
 
         %Deal with boundary conditions:
         
-        disp('      boundary conditions for diagnostic model: ');
-        %Build gridonicefront, array of boundary grids belonging to the icefront:
-        gridinsideicefront=ArgusContourToMesh(md.elements,md.x,md.y,expread('Front.exp',1),'node',2);
-        gridonicefront=double(md.gridonboundary & gridinsideicefront);
+        disp('      boundary conditions');
+        md=SetIceShelfBC(md,'Front.exp');
 
-        md.gridondirichlet_diag=zeros(md.numberofgrids,1);
-        pos=find(md.gridonboundary & ~gridonicefront);md.gridondirichlet_diag(pos)=1;
-        md.dirichletvalues_diag=100*ones(md.numberofgrids,2)*[0,0;0,1];
-
-        pos=find(gridonicefront(md.segments(:,1)) | gridonicefront(md.segments(:,2)));
-        md.segmentonneumann_diag=md.segments(pos,:);
-        md.neumannvalues_diag=NaN*ones(length(md.segmentonneumann_diag),1); %dynamic boundary conditions (water pressure)
-
-        disp('      boundary conditions for prognostic model: ');
-        md.gridondirichlet_prog=zeros(md.numberofgrids,1);
-        md.dirichletvalues_prog=zeros(md.numberofgrids,1);
-        pos=find(gridonicefront(md.segments(:,1)) | gridonicefront(md.segments(:,2)));
-        md.segmentonneumann_prog=md.segments(pos,:);
-        md.neumannvalues_prog=zeros(size(md.segmentonneumann_prog,1),1);
-        md.neumannvalues_prog(:)=NaN; %free radiation
-        
-        pos=find(gridonicefront(md.segments(:,1)) | gridonicefront(md.segments(:,2)));
-        md.segmentonneumann_prog2=md.segments(pos,:);
-        md.neumannvalues_prog2=zeros(size(md.segmentonneumann_prog2,1),1);
-        md.neumannvalues_prog2(:)=NaN; %free radiation
-        
-        disp('      boundary conditions for thermal model: ');
-        md.gridondirichlet_thermal=ones(md.numberofgrids,1); %surface temperature
-        md.dirichletvalues_thermal=md.observed_temperature;
-        md.geothermalflux=zeros(md.numberofgrids,1); 
-
-
-% Some Cielo code, ignore.
-if strcmp(md.cluster,'yes')
-        ServerDisconnect;
-end   

issm/trunk/examples/Bumps/Bump1_surface_bed/runme.m

@@ -1 @@
-% This file can be run to see the response to a bump on the bed of an ice sheet. The bump is not 
-% transmitted to the surface of the ice. The results of the three models are saved in order to be 
-% compared. The geometry is square and the ice sheet is flat. The velocity of the icesheet is 
-% 100 m/yr in one particular direction on every border. Just run this file in Matlab, i
-% with a properly setup Ice code. 
-% Take care to have the same density and number of layers to compare the results easily.
-
-% Create macayeal model;
-
-mdm=model;
-mdm=mesh(mdm,'DomainOutline.exp',150000);
-mdm=geography(mdm,'','');
-mdm=parameterize(mdm,'Square.par');
-mdm=extrude(mdm,6,3);
-mdm=setelementstype(mdm,'macayeal','all');
-
-%Compute solution with Ice model
-mdm=solve(mdm,'diagnostic','ice');
-
-save modelmacayeal mdm
-
-% Create pattyn model;
-
-mdp=model;
-mdp=mesh(mdp,'DomainOutline.exp',150000);
-mdp=geography(mdp,'','');
-mdp=parameterize(mdp,'Square.par');
-mdp=extrude(mdp,6,3);
-mdp=setelementstype(mdp,'pattyn','all');
-
-%Compute solution with Ice model
-mdp=solve(mdp,'diagnostic','ice');
-
-save modelpattyn mdp
-
-% Create stokes model;
-
-mds=model;
-mds=mesh(mds,'DomainOutline.exp',150000);
-mds=geography(mds,'','');
-mds=parameterize(mds,'Square.par');
-mds=extrude(mds,6,3);
-mds=setelementstype(mds,'pattyn','all','stokes','all');
-
-%Compute solution with Ice model
-mds=solve(mds,'diagnostic','ice');
-
-save modelstokes mds
+if 0,
+        md=mesh(model,'DomainOutline.exp',15000);
+        md=geography(md,'','');
+        md=parameterize(md,'Square.par');
+        md=extrude(md,2,3);
+        md=setelementstype(md,'stokes','all');
+        md.cluster.np=14;
+        md.verbose.convergence=1;
+        md.spcvelocity(:,5)=1000;
+        md=solve(md,'analysis_type',DiagnosticSolutionEnum);
+else
+        md.vx=PatchToVec(md.results.DiagnosticSolution.Vx);
+        md.vy=PatchToVec(md.results.DiagnosticSolution.Vy);
+        md.accumulation_rate(:)=0;
+        md.melting_rate(:)=0;
+        md.thickness(:)=1;
+        md.bed=md.surface-md.thickness;
+        md.dt=1;
+        md=solve(md,'analysis_type',PrognosticSolutionEnum);
+        beep;pause(0.2);beep;
+        plotmodel(md,'data',PatchToVec(md.results.PrognosticSolution.Thickness)-md.thickness) 
+end

issm/trunk/externalpackages/export_fig/ghostscript.m

@@ -124 +124 @@
 
 function current_gs_path_str = current_gs_path
-current_gs_path_str = 'C:\Program Files\gs\gs8.71\bin\gswin32c.exe';
+current_gs_path_str = 'gs';
 return

issm/trunk/src/m/utils/Mesh/FixMesh.m

@@ -15 +15 @@
 y2=y;
 value2=value;
-
 
 %First, look for orphan vertices, and take them out.
@@ -39 +38 @@
 end
 
-
 %Check all triangles are well oriented.
 aires=GetAreas(index2,x2,y2);

issm/trunk/src/mex/Test/Test.cpp

@@ -11 +11 @@
         char    *name    = NULL;
         bool     logical;
+        double  *matrix    = NULL;
+        int      numel;
 
         printf("Parsing options:\n");
@@ -24 +26 @@
         printf("\n value of \"CELL[0,2]\" is %g\n",test);
 
-        options->Get(&logical,"logical");
+        options->Get(&logical,"logical",true);
         printf("\n value of \"logical\" is %s\n",logical?"true":"false");
+
+        options->Get(&matrix,&numel,"matrix");
+        printarray(matrix,numel);
 
         delete options;
         xfree((void**)&name);
+        xfree((void**)&matrix);
 }