Changeset 11191


Ignore:
Timestamp:
01/24/12 14:01:23 (13 years ago)
Author:
seroussi
Message:

finally finished to fix the NR !

Location:
issm/trunk-jpl/test
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/test/NightlyRun/test1108.m

    r11031 r11191  
    66
    77for i=1:length(L_list),
    8         L=3*L_list{i};
    9         nx=60; %numberof nodes in x direction
    10         ny=60;
     8        L=L_list{i};
     9        nx=30; %numberof nodes in x direction
     10        ny=30;
    1111        md=model;
    1212        md=squaremesh(md,L,L,nx,ny);
    1313        md=setmask(md,'',''); %ice sheet test
    1414        md=parameterize(md,'../Par/ISMIPD.par');
    15         md=extrude(md,6,1);
     15        md=extrude(md,10,1);
    1616
    17         md=setflowequation(md,'stokes','all');
     17        md=setflowequation(md,'pattyn','all');
    1818
    1919        %We need one grd on dirichlet: the 4 corners are set to zero
    20         %md.diagnostic.spcvx=NaN*ones(md.mesh.numberofvertices,1);
    21         %md.diagnostic.spcvy=NaN*ones(md.mesh.numberofvertices,1);
    22         %md.diagnostic.spcvz=NaN*ones(md.mesh.numberofvertices,1);
     20        md.diagnostic.spcvx=NaN*ones(md.mesh.numberofvertices,1);
     21        md.diagnostic.spcvy=NaN*ones(md.mesh.numberofvertices,1);
     22        md.diagnostic.spcvz=NaN*ones(md.mesh.numberofvertices,1);
    2323       
    24         %pos=find((md.mesh.x==0 | md.mesh.x==max(md.mesh.x)) & (md.mesh.y==0 | md.mesh.y==max(md.mesh.y)));
    25         %md.diagnostic.spcvx(pos)=0;
    26         %md.diagnostic.spcvy(pos)=0;
    27         %md.diagnostic.spcvz(pos)=0;
     24        pos=find(md.mesh.vertexonbed & (md.mesh.x==0 | md.mesh.x==max(md.mesh.x)) & (md.mesh.y==0 | md.mesh.y==max(md.mesh.y)));
     25        md.diagnostic.spcvx(pos)=0;
     26        md.diagnostic.spcvy(pos)=0;
     27        md.diagnostic.spcvz(pos)=0;
    2828
    29         %%Create MPCs to have periodic boundary conditions
    30         %posx=find(md.mesh.x==0);
    31         %posx2=find(md.mesh.x==max(md.mesh.x));
     29        %Create MPCs to have periodic boundary conditions
     30        posx=find(md.mesh.x==0);
     31        posx2=find(md.mesh.x==max(md.mesh.x));
    3232
    33         %posy=find(md.mesh.y==0 & md.mesh.x~=0 & md.mesh.x~=max(md.mesh.x)); %Don't take the same nodes two times
    34         %posy2=find(md.mesh.y==max(md.mesh.y) & md.mesh.x~=0 & md.mesh.x~=max(md.mesh.x));
     33        posy=find(md.mesh.y==0 & md.mesh.x~=0 & md.mesh.x~=max(md.mesh.x)); %Don't take the same nodes two times
     34        posy2=find(md.mesh.y==max(md.mesh.y) & md.mesh.x~=0 & md.mesh.x~=max(md.mesh.x));
    3535
    36         %md.diagnostic.vertex_pairing=[posx,posx2;posy,posy2];
     36        md.diagnostic.vertex_pairing=[posx,posx2;posy,posy2];
    3737
    3838        %Compute the diagnostic
    3939        md.cluster=generic('name',oshostname(),'np',8);
     40        md.verbose=verbose('convergence',true);
     41        md=solve(md,DiagnosticSolutionEnum);
     42        md.diagnostic.reltol=NaN;
     43        md.diagnostic.abstol=NaN;
     44        md.diagnostic.vertex_pairing=[];
     45        %We need one grd on dirichlet: the 4 corners are set to zero
     46        md.diagnostic.spcvx=NaN*ones(md.mesh.numberofvertices,1);
     47        md.diagnostic.spcvy=NaN*ones(md.mesh.numberofvertices,1);
     48        md.diagnostic.spcvz=NaN*ones(md.mesh.numberofvertices,1);
     49        pos=find(md.mesh.y==0 | md.mesh.x==0 | md.mesh.x==max(md.mesh.x) | md.mesh.y==max(md.mesh.y)); %Don't take the same nodes two times
     50        md.diagnostic.spcvx(pos)=md.results.DiagnosticSolution.Vx(pos);
     51        md.diagnostic.spcvy(pos)=md.results.DiagnosticSolution.Vy(pos);
     52        md=setflowequation(md,'stokes','all');
    4053        md=solve(md,DiagnosticSolutionEnum);
    4154
     
    5770        'Vx80km','Vy80km','Vz80km',...
    5871        'Vx160km','Vy160km','Vz160km'
    59 }
     72};
    6073field_tolerances={...
    6174        1e-11,1e-11,1e-11,...
Note: See TracChangeset for help on using the changeset viewer.