Changeset 4869


Ignore:
Timestamp:
07/29/10 09:29:57 (15 years ago)
Author:
seroussi
Message:

eismint tests (shelves)

Location:
issm/trunk/test/Validation/EISMINT
Files:
9 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk/test/Validation/EISMINT/Diagnostic1/Square.par

    r3768 r4869  
    66md.thickness=500*ones(md.numberofgrids,1);
    77md.firn_layer=0*ones(md.numberofgrids,1);
    8 md.bed=-di*md.thickness;
     8md.bed=-md.rho_ice/md.rho_water*md.thickness;
    99md.surface=md.bed+md.thickness;
    1010
    … …  
    1515pos=find(md.elementoniceshelf);
    1616md.drag_coefficient(md.elements(pos,:))=0;
    17 md.drag-p=ones(md.numberofelements,1);
     17md.drag_p=ones(md.numberofelements,1);
    1818md.drag_q=ones(md.numberofelements,1);
    1919

  • issm/trunk/test/Validation/EISMINT/Diagnostic1/runme.m

    r899 r4869  
    88md=parameterize(md,'Square.par');
    99md=setelementstype(md,'macayeal','all'); %MacAyeal's model and 2d
     10md.verbose=10;
    1011
    1112%Compute solution for MacAyeal's model
    12 md=solve(md,'analysis_type','diagnostic','package','ice');
     13md=solve(md,'analysis_type',DiagnosticSolutionEnum);
    1314
    1415%plot results
    15 plotmodel(md,'data','vx','data','vy','contourlevels',{0,20,40,-20,-40},'contourlevels',{-100,-200,-300,-400,-500,-600})
     16vx=zeros(md.numberofgrids,1);
     17vx(md.results.DiagnosticSolution.Vx.index)=md.results.DiagnosticSolution.Vx.value;
     18vy=zeros(md.numberofgrids,1);
     19vy(md.results.DiagnosticSolution.Vy.index)=md.results.DiagnosticSolution.Vy.value;
     20
     21plotmodel(md,'data',vx,'data',vy,'contourlevels',{0,20,40,-20,-40},'contourlevels',{-100,-200,-300,-400,-500,-600})

  • issm/trunk/test/Validation/EISMINT/Diagnostic2/Square.par

    r3768 r4869  
    66md.thickness=500*ones(md.numberofgrids,1);
    77md.firn_layer=0*ones(md.numberofgrids,1);
    8 md.bed=-di*md.thickness;
     8md.bed=-md.rho_ice/md.rho_water*md.thickness;
    99md.surface=md.bed+md.thickness;
    1010
    … …  
    3232md=SetMarineIceSheetBC(md,'Front.exp');
    3333pos=find(md.y==max(md.y)); %grids on the upper boundary condition
    34 md.dirichletvalues_diag(pos,2)=400*(((md.x(pos)-100000)/25000).^2-ones(size(pos,1),1)).*heaviside((1+eps)*ones(size(pos,1),1)-((md.x(pos)-100000)/25000).^2);
     34md.spcvelocity(pos,5)=400*(((md.x(pos)-100000)/25000).^2-ones(size(pos,1),1)).*heaviside((1+eps)*ones(size(pos,1),1)-((md.x(pos)-100000)/25000).^2);
     35
     36disp('      cluster properties');
     37md.np=3;
     38md.cluster='wilkes';

  • issm/trunk/test/Validation/EISMINT/Diagnostic2/runme.m

    r899 r4869  
    88md=parameterize(md,'Square.par');
    99md=setelementstype(md,'macayeal','all'); %MacAyeal's model and 2d
     10md.verbose=10;
    1011
    1112%Compute solution for MacAyeal's model
    12 md=solve(md,'analysis_type','diagnostic','package','ice');
     13md=solve(md,'analysis_type',DiagnosticSolutionEnum);
    1314
    1415%plot results
    15 plotmodel(md,'data','vx','data','vy','contourlevels',{0,20,40,-20,-40,-60,60},'contourlevels',{-100,-200,-300,-400,-500,-600,-700,-800})
     16vx=zeros(md.numberofgrids,1);
     17vx(md.results.DiagnosticSolution.Vx.index)=md.results.DiagnosticSolution.Vu.value;
     18vy=zeros(md.numberofgrids,1);
     19vy(md.results.DiagnosticSolution.Vy.index)=md.results.DiagnosticSolution.Vy.value;
     20plotmodel(md,'data',vx,'data',vy,'contourlevels',{0,20,40,-20,-40,-60,60},'contourlevels',{-100,-200,-300,-400,-500,-600,-700,-800})

  • issm/trunk/test/Validation/EISMINT/MassConservation/Square.par

    r3768 r4869  
    1515pos=find(md.elementoniceshelf);
    1616md.drag_coefficient(md.elements(pos,:))=0;
    17 md.drag-p=ones(md.numberofelements,1);
     17md.drag_p=ones(md.numberofelements,1);
    1818md.drag_q=ones(md.numberofelements,1);
    1919
    … …  
    4848md.spcthickness(pos,1)=1;
    4949md.spcthickness(pos,2)=500;
     50md.prognostic_DG=1;

  • issm/trunk/test/Validation/EISMINT/MassConservation/runme.m

    r3573 r4869  
    2828        i=i+1;
    2929        time=time+md.dt;
    30         md=solve(md,'analysis_type','prognostic2','package','ice');
    31         %thickness=full(sparse(reshape(md.elements',3*md.numberofelements,1),1,md.results.prognostic2.thickness));
    32         %md.thickness=thickness./connectivity;
    33         md.thickness=md.results.prognostic2.thickness;
     30        md=solve(md,'analysis_type',PrognosticSolutionEnum);
     31        thickness   =zeros(md.numberofgrids,1);
     32        connectivity=sparse(md.results.PrognosticSolution.Thickness.index(:),1,1);
     33        thickness   =sparse(md.results.PrognosticSolution.Thickness.index(:),1,md.results.PrognosticSolution.Thickness.value(:));
     34        md.thickness=full(thickness./connectivity);
    3435        md.thickness(pos)=500+100*sin(2*pi*time/200)*ones(size(pos,1),1);
    3536        md.surface=md.bed+md.thickness;

  • issm/trunk/test/Validation/EISMINT/MassConservation/runme2.m

    r3581 r4869  
    3333        i=i+1;
    3434        time=time+md.dt;
    35         md=solve(md,'analysis_type','prognostic2','package','ice');
    36         %thickness=full(sparse(reshape(md.elements',3*md.numberofelements,1),1,md.results.prognostic2.thickness));
    37         %md.thickness=thickness./connectivity;
    38         md.thickness=md.results.prognostic2.thickness;
     35        md=solve(md,'analysis_type',PrognosticSolutionEnum);
     36        thickness   =zeros(md.numberofgrids,1);
     37        connectivity=sparse(md.results.PrognosticSolution.Thickness.index(:),1,1);
     38        thickness   =sparse(md.results.PrognosticSolution.Thickness.index(:),1,md.results.PrognosticSolution.Thickness.value(:));
     39        md.thickness=full(thickness./connectivity);
    3940        md.thickness(posLEFT)=500+100*sin((2*pi/200)*(time-sinus*md.y(posLEFT)/sqrt(md.vx(1)^2+md.vy(1)^2)));
    4041        md.thickness(posDOWN)=500+100*sin((2*pi/200)*(time-cosinus*md.x(posDOWN)/sqrt(md.vx(1)^2+md.vy(1)^2)));

  • issm/trunk/test/Validation/EISMINT/Transient/Square.par

    r3768 r4869  
    66md.thickness=500*ones(md.numberofgrids,1);
    77md.firn_layer=0*ones(md.numberofgrids,1);
    8 md.bed=-di*md.thickness;
     8md.bed=-md.rho_ice/md.rho_water*md.thickness;
    99md.surface=md.bed+md.thickness;
    1010
    … …  
    2626
    2727disp('      creating accumulation rates');
    28 md.accumulation_rate=0.2*ones(md.numberofgrids,1); %0.2m/a
     28md.accumulation_rate=200*ones(md.numberofgrids,1); %0.2m/a
    2929md.melting_rate=0*ones(md.numberofgrids,1); %0m/a
    3030
    3131disp('      boundary conditions ');
    3232md=SetMarineIceSheetBC(md,'Front.exp');
     33pos=find(md.y==max(md.y)); %grids on the upper boundary condition
     34md.spcvelocity(pos,5)=400*(((md.x(pos)-100000)/25000).^2-ones(size(pos,1),1)).*heaviside((1+eps)*ones(size(pos,1),1)-((md.x(pos)-100000)/25000).^2);
    3335
    34 disp('      boundary conditions for prognostic model: ');
    35 md.gridondirichlet_prog=zeros(md.numberofgrids,1);
    36 md.dirichletvalues_prog=zeros(md.numberofgrids,1);
     36disp('      cluster properties');
     37md.np=3;
     38md.cluster='wilkes';
    3739
    38 pos=find(md.y==max(md.y)); %grids on the upper boundary condition
    39 md.dirichletvalues_diag(pos,2)=400*(((md.x(pos)-100000)/25000).^2-ones(size(pos,1),1)).*heaviside((1+eps)*ones(size(pos,1),1)-((md.x(pos)-100000)/25000).^2);
    40 md.neumannvalues_prog=zeros(size(md.segmentonneumann_prog,1),1);
    41 md.neumannvalues_prog(:)=NaN; %free radiation
    42 
    43 pos=find(md.y==max(md.y)); %grids on the upper boundary condition
    44 md.gridondirichlet_prog(pos)=1;
    45 md.dirichletvalues_prog(pos)=500;
     40disp('      prognostic boundary conditions ');
     41pos=find(md.y==max(md.y));
     42md.spcthickness(pos,1)=1;
     43md.spcthickness(pos,2)=md.thickness(pos);

  • issm/trunk/test/Validation/EISMINT/Transient/runme.m

    r899 r4869  
    33%Evolution of the ice shelf
    44md=model;
    5 md=mesh(md,'DomainOutline.exp',7100);
     5md=mesh(md,'DomainOutline.exp',5100);
    66md=geography(md,'all','');
    77md=parameterize(md,'Square.par');
    88md=setelementstype(md,'macayeal','all'); %MacAyeal's model and 2d
    99
     10md.np=5;
     11md.cluster='wilkes';
     12md.verbose=3;
    1013%Compute solution for MacAyeal's model
    11 md=solve(md,'analysis_type','diagnostic','package','ice');
     14md=solve(md,'analysis_type',DiagnosticSolutionEnum);
    1215
    1316%Launch transient solution
    14 md.ndt=5000*md.yts; %time of evolution
    15 md.dt=500*md.yts; %time vetween two iterations
    16 md=solve(md,'analysis_type','transient');
     17vx=zeros(md.numberofgrids,1);
     18vx(md.results.DiagnosticSolution.Vx.index)=md.results.DiagnosticSolution.Vx.value;
     19vy=zeros(md.numberofgrids,1);
     20vy(md.results.DiagnosticSolution.Vy.index)=md.results.DiagnosticSolution.Vy.value;
     21md.vx=vx;
     22md.vy=vy;
     23
     24md.ndt=1000*md.yts; %time of evolution
     25md.dt=20*md.yts; %time vetween two iterations
     26md=solve(md,'analysis_type',Transient2DSolutionEnum);
    1727
    1828%plot results
    19 plotmodel(md,'data',md.results.transient(end).vx,'data',md.results.transient(end).vy,'data',md.results.transient(end).thickness)
     29plotmodel(md,'data',md.results.Transient2DSolution(end).Vx,'data',md.results.Transient2DSolution(end).Vy,'data',md.results.Transient2DSolution(end).Thickness)
Note: See TracChangeset for help on using the changeset viewer.