Index: /issm/trunk/src/m/classes/public/BuildQueueingScript.m =================================================================== --- /issm/trunk/src/m/classes/public/BuildQueueingScript.m (revision 1110) +++ /issm/trunk/src/m/classes/public/BuildQueueingScript.m (revision 1111) @@ -9,5 +9,13 @@ %First try and figure out if there is a special script for thie particular cluster function_name=['BuildQueueingScript' md.cluster]; + +%some specific treatment of identical cluster, gemini, castor and pollux +if strcmpi(md.cluster,'castor') || strcmpi(md.cluster,'pollux'), + function_name='BuildQueueingScriptgemini'; +end + if exist(function_name,'file'), + + %Call this function: eval([function_name '(md,executionpath,codepath);']); Index: /issm/trunk/src/m/classes/public/BuildQueueingScriptGeneric.m =================================================================== --- /issm/trunk/src/m/classes/public/BuildQueueingScriptGeneric.m (revision 1110) +++ /issm/trunk/src/m/classes/public/BuildQueueingScriptGeneric.m (revision 1111) @@ -15,20 +15,4 @@ fprintf(fid,'#!/bin/sh\n'); fprintf(fid,'rm -rf %s/%s.lock\n',executionpath,md.name); -fprintf(fid,'mpirun -np %i ',md.np); - -if strcmpi(md.analysis_type,'diagnostic'), - fprintf(fid,'%s/diagnostic.exe',codepath); -elseif strcmpi(md.analysis_type,'control'), - fprintf(fid,'%s/control.exe',codepath); -elseif strcmpi(md.analysis_type,'thermal'), - fprintf(fid,'%s/thermal.exe',codepath); -elseif strcmpi(md.analysis_type,'prognostic'), - fprintf(fid,'%s/prognostic.exe',codepath); -elseif strcmpi(md.analysis_type,'transient'), - fprintf(fid,'%s/transient.exe',codepath); -else - error('BuildQueueingScriptGeneric error message: unsupported solution type!'); -end - -fprintf(fid,' %s %s.bin %s.outbin %s.lock 2> %s.errlog >%s.outlog & ',executionpath,md.name,md.name,md.name,md.name,md.name); +fprintf(fid,'mpirun -np %i %s/%s.exe %s %s.bin %s.outbin %s.lock 2> %s.errlog >%s.outlog & ',md.np,codepath,md.analysis_type,executionpath,md.name,md.name,md.name,md.name,md.name); fclose(fid); Index: /issm/trunk/src/m/classes/public/BuildQueueingScriptgemini.m =================================================================== --- /issm/trunk/src/m/classes/public/BuildQueueingScriptgemini.m (revision 1110) +++ /issm/trunk/src/m/classes/public/BuildQueueingScriptgemini.m (revision 1111) @@ -1,3 +1,3 @@ -function BuildQueueingScriptgemini(md,solutiontype,executionpath,codepath) +function BuildQueueingScriptgemini(md,executionpath,codepath) %BUILDQUEUEINGSCRIPTGEMINI - ... % @@ -12,22 +12,17 @@ end -fprintf(fid,'#!/bin/sh\n'); -fprintf(fid,'#BSUB -n %i -W %i\n',md.np,md.time); -fprintf(fid,'#BSUB -J %s\n',md.name); +fprintf(fid,'#!/bin/csh\n'); +fprintf(fid,'#PBS -l walltime=%i\n',md.time); +fprintf(fid,'#PBS -l ncpus=%i\n',md.np); if ~isempty(md.queue), - fprintf(fid,'#BSUB -q %s\n',md.queue); + fprintf(fid,'#PBS -q %s\n',md.queue); end -fprintf(fid,'#BSUB -o %s.outlog -e %s.errlog\n',md.name,md.name); -fprintf(fid,'cd %s\n',executionpath); -fprintf(fid,'rm -rf %s.outlog %s.errlog %s.lock\n',md.name,md.name,md.name); +fprintf(fid,'#PBS -o %s.outlog \n',md.name); +fprintf(fid,'#PBS -e %s.errlog \n',md.name); -if strcmpi(solutiontype,'diagnostic') , - fprintf(fid,'mpirun -np %i %s/diagnostic.exe %s %s.bin %s.outbin %s.lock',md.np,codepath,executionpath,md.name,md.name,md.name); -elseif strcmpi(solutiontype,'control') , - fprintf(fid,'mpirun -np %i %s/control.exe %s %s.bin %s.outbin %s.lock',md.np,codepath,executionpath,md.name,md.name,md.name); -elseif strcmpi(solutiontype,'thermal') , - fprintf(fid,'mpirun -np %i %s/thermal.exe %s %s.bin %s.outbin %s.lock',md.np,codepath,executionpath,md.name,md.name,md.name); -else - error('BuildQueueingScriptgemini error message: unsupported solution type!'); -end +fprintf(fid,'setenv PBS_O_WORKDIR %s\n',executionpath); +fprintf(fid,'cd $PBS_O_WORKDIR\n'); +fprintf(fid,'setenv OMP_NUM_THREADS 1\n'); +fprintf(fid,'mpirun -np %i %s/%s.exe %s %s.bin %s.outbin %s.lock',md.np,codepath,md.analysis_type,executionpath,md.name,md.name,md.name); + fclose(fid); Index: /issm/trunk/src/m/classes/public/LaunchQueueJob.m =================================================================== --- /issm/trunk/src/m/classes/public/LaunchQueueJob.m (revision 1110) +++ /issm/trunk/src/m/classes/public/LaunchQueueJob.m (revision 1111) @@ -7,4 +7,12 @@ %First try and figure out if there is a special script for thie particular cluster function_name=['LaunchQueueJob' md.cluster]; + +%some specific treatment of identical cluster, gemini, castor and pollux +if strcmpi(md.cluster,'castor') || strcmpi(md.cluster,'pollux'), + function_name='LaunchQueueJobgemini'; +end + + + if exist(function_name,'file'), %Call this function: Index: /issm/trunk/src/m/classes/public/LaunchQueueJobgemini.m =================================================================== --- /issm/trunk/src/m/classes/public/LaunchQueueJobgemini.m (revision 1110) +++ /issm/trunk/src/m/classes/public/LaunchQueueJobgemini.m (revision 1111) @@ -16,7 +16,17 @@ %jpload both files to cluster disp('uploading input file, queueing script and variables script'); -eval(['!scp ' md.name '.bin' ' ' md.name '.queue ' md.cluster ':' executionpath]); +if md.qmu_analysis, + eval(['!scp ' md.name '.bin' ' ' md.name '.queue ' md.name '.qmu.in ' md.cluster ':' executionpath]); +else + eval(['!scp ' md.name '.bin' ' ' md.name '.queue ' md.cluster ':' executionpath]); +end + +%new gemini cannot launch across cluster using ssh +disp(['launch solution sequence on remote cluster by logging into it and typing qsub < ' md.name '.queue']); +return; disp('launching solution sequence on remote cluster'); + %now call the queuing script to launch the job. -system(['ssh ' md.cluster ' ''cd ' executionpath ' && bsub < ' md.name '.queue ''']); +system(['ssh ' md.cluster ' ''cd ' executionpath ' && qsub < ' md.name '.queue ''']); + Index: /issm/trunk/src/m/classes/public/ThicknessCorrection.m =================================================================== --- /issm/trunk/src/m/classes/public/ThicknessCorrection.m (revision 1110) +++ /issm/trunk/src/m/classes/public/ThicknessCorrection.m (revision 1111) @@ -34,4 +34,10 @@ %modify thickness for i=1:length(pos_shelf) + + if mod(i,100)==0, + fprintf('\b\b\b\b\b\b\b') + fprintf('%5.2f%s',i/length(pos_shelf)*100,' %'); + end + %search the grid on ice sheet the closest to i [d posd]=min(sqrt((md.x(pos_shelf(i))-md.x(pos_GL)).^2+(md.y(pos_shelf(i))-md.y(pos_GL)).^2)); Index: /issm/trunk/src/m/classes/public/displayqmu.m =================================================================== --- /issm/trunk/src/m/classes/public/displayqmu.m (revision 1110) +++ /issm/trunk/src/m/classes/public/displayqmu.m (revision 1111) @@ -89,3 +89,5 @@ if isempty(md.dakotadat), disp(sprintf(' dakotadat: N/A')); else disp(sprintf(' dakotadat: [%ix%i] (can be accessed by typing md.dakotadat)',size(md.dakotadat)));end disp(sprintf(' npart : %i (number of partitions for semi-descrete qmu)',md.npart)); + disp(sprintf(' part : [%i] (user provided mesh partition, defaults to metis if not specified)',length(md.part))); + end Index: /issm/trunk/src/m/classes/public/geography2.m =================================================================== --- /issm/trunk/src/m/classes/public/geography2.m (revision 1111) +++ /issm/trunk/src/m/classes/public/geography2.m (revision 1111) @@ -0,0 +1,115 @@ +function md=geography2(md,landname,iceshelfname,icesheetname) + +%Get assigned fields +x=md.x; +y=md.y; +elements=md.elements; + +%recover elements and grids on land. +[gridonland,elementonland]=ContourToMesh(elements,x,y,expread(landname,1),'element and node',1); + +%use segments to make sure gridonland is correct at the boundary: +gridonland(md.segments(:,1))=1; +gridonland(md.segments(:,2))=1; + +%now, any elements with 3 grids on land should be on land: +elementsonwater=find(~elementonland); +wrongelements=elementsonwater(find(( gridonland(md.elements(elementsonwater,1)) + gridonland(md.elements(elementsonwater,2)) + gridonland(md.elements(elementsonwater,3)) ... + )==3)); +elementonland(wrongelements)=1; + +%finally, any element with its barycentre in the water should be in the water! +xelem=x(md.elements)*[1;1;1]/3; yelem=y(md.elements)*[1;1;1]/3; +baryonland=ContourToNodes(xelem,yelem,expread(landname,1),1); +pos=find(~baryonland); elementonland(pos)=0; + +%figure out which elements on land are actually in the middle of the ocean! +pos1=find(elementonland); +connectedtoland=md.elementconnectivity(pos1,:); +pos=find(connectedtoland); connectedtoland(pos)=1-elementonland(connectedtoland(pos)); +connectedtolandsum=sum(connectedtoland,2); +waterelements=pos1(find(connectedtolandsum==3)); + +%now use waterelements to correct elementonland: +elementonland(waterelements)=0; + +%recover arrays of ice shelf grids and elements, and ice sheet grids and elements. +if strcmp(iceshelfname,''), %no iceshelf contour file, we are dealing with a pure ice sheet. + gridoniceshelf=zeros(md.numberofgrids,1); + elementoniceshelf=zeros(md.numberofelements,1); +elseif strcmp(iceshelfname,'all'), %we are dealing with a pure ice shelf. + gridoniceshelf=ones(md.numberofgrids,1); + elementoniceshelf=ones(md.numberofelements,1); +else + [gridoniceshelf,elementoniceshelf]=ContourToMesh(elements,x,y,expread(iceshelfname,1),'element and node',1); +end + +if strcmp(icesheetname,''), %no icesheet contour file, we are dealing with a pure ice shelf. + gridonicesheet=zeros(md.numberofgrids,1); + elementonicesheet=zeros(md.numberofelements,1); +else + [gridonicesheet,elementonicesheet]=ContourToMesh(elements,x,y,expread(icesheetname,1),'element and node',1); +end + + +%now correct, so that none of the iceshelf and icesheet elements and nodes are in the water. +pos=find(~elementonland); +elementoniceshelf(pos)=0; +elementonicesheet(pos)=0; + +pos=find(~gridonland); +gridoniceshelf(pos)=0; +gridonicesheet(pos)=0; + +%create gridonwater and elementonwater: +gridonwater=double(~gridonland); +elementonwater=double(~elementonland); + +%correct for islands: +gridoniceshelf=double(gridoniceshelf & ~gridonicesheet); +elementoniceshelf=double(elementoniceshelf & ~elementonicesheet); + +%now, icesheets are everything except iceshelves and water +gridonicesheet=double(~gridoniceshelf & ~gridonwater); +elementonicesheet=double(~elementoniceshelf & ~elementonwater); + +%Return: +md.gridoniceshelf=gridoniceshelf; +md.elementoniceshelf=elementoniceshelf; + +md.gridonwater=gridonwater; +md.elementonwater=elementonwater; + +md.gridonicesheet=gridonicesheet; +md.elementonicesheet=elementonicesheet; + +%Deal with segments on neumann: +pos=find(md.segmentmarkers==2); %internal segments. +segments=md.segments(pos,:); + +%swap segments so they are all oriented the same way +for i=1:length(segments), + pair=find( (md.elements(:,1)==segments(i,1) | md.elements(:,2)==segments(i,1) | md.elements(:,3)==segments(i,1) ) & ... + (md.elements(:,1)==segments(i,2) | md.elements(:,2)==segments(i,2) | md.elements(:,3)==segments(i,2) ) ... + ); + if md.elementonwater(pair(1)), + segments(i,3)=pair(2); + else + segments(i,3)=pair(1); + end + ord1=find(segments(i,1)==md.elements(segments(i,3),:)); + ord2=find(segments(i,2)==md.elements(segments(i,3),:)); + + %swap segment grids if necessary + if ( (ord1==1 & ord2==2) | (ord1==2 & ord2==3) | (ord1==3 & ord2==1) ), + temp=segments(i,1); + segments(i,1)=segments(i,2); + segments(i,2)=temp; + end + segments(i,1:2)=fliplr(segments(i,1:2)); +end +md.segmentonneumann_diag=segments; +md.counter=2; +md.segmentmarkers(:)=1; + + Index: /issm/trunk/src/m/classes/public/ismodelselfconsistent.m =================================================================== --- /issm/trunk/src/m/classes/public/ismodelselfconsistent.m (revision 1110) +++ /issm/trunk/src/m/classes/public/ismodelselfconsistent.m (revision 1111) @@ -162,4 +162,19 @@ disp(['concurrency should be set to 1 when running dakota in library mode']); bool=0;return; + end + if ~isempty(md.part), + if numel(md.part)~=md.numberofgrids, + disp(['user supplied partition for qmu analysis should have size md.numberofgrids x 1 ']); + bool=0;return; + end + if find(md.part)>=md.numberofgrids, + disp(['user supplied partition should be indexed from 0 (c-convention)']); + bool=0;return; + end + if md.npart~=md.numberofgrids, + disp(['user supplied partition should have same size as md.npart']); + bool=0;return; + end + end end Index: /issm/trunk/src/m/classes/public/loadresultsfromcluster.m =================================================================== --- /issm/trunk/src/m/classes/public/loadresultsfromcluster.m (revision 1110) +++ /issm/trunk/src/m/classes/public/loadresultsfromcluster.m (revision 1111) @@ -24,7 +24,7 @@ else if md.qmu_analysis, - system(['scp ' md.cluster ':' executionpath '/' md.name '.outlog ' executionpath '/' md.name '.errlog ' executionpath '/' md.name '.outbin ' executionpath '/' md.name '.qmu.err ' executionpath '/' md.name '.qmu.out ' executionpath '/dakota_tabular.dat ./']); %get outlog, errlog and outbin files + system(['scp ' md.cluster ':' executionpath '/' md.name '.outlog ' md.cluster ':' executionpath '/' md.name '.errlog ' md.cluster ':' executionpath '/' md.name '.outbin ' md.cluster ':' executionpath '/' md.name '.qmu.err ' md.cluster ':' executionpath '/' md.name '.qmu.out ' md.cluster ':' executionpath '/dakota_tabular.dat ./']); %get outlog, errlog and outbin files else - system(['scp ' md.cluster ':' executionpath '/' md.name '.outlog ' executionpath '/' md.name '.errlog ' executionpath '/' md.name '.outbin ./']); %get outlog, errlog and outbin files + system(['scp ' md.cluster ':' executionpath '/' md.name '.outlog ' md.cluster ':' executionpath '/' md.name '.errlog ' md.cluster ':' executionpath '/' md.name '.outbin ./']); %get outlog, errlog and outbin files end end Index: /issm/trunk/src/m/classes/public/marshall.m =================================================================== --- /issm/trunk/src/m/classes/public/marshall.m (revision 1110) +++ /issm/trunk/src/m/classes/public/marshall.m (revision 1111) @@ -73,4 +73,5 @@ WriteData(fid,md.elementoniceshelf,'Mat','elementoniceshelf'); +WriteData(fid,md.elementonwater,'Mat','elementonwater'); WriteData(fid,md.gridonicesheet,'Mat','gridonicesheet'); WriteData(fid,md.gridoniceshelf,'Mat','gridoniceshelf'); Index: /issm/trunk/src/m/classes/public/mesh.m =================================================================== --- /issm/trunk/src/m/classes/public/mesh.m (revision 1110) +++ /issm/trunk/src/m/classes/public/mesh.m (revision 1111) @@ -77,4 +77,8 @@ md.elementonsurface=ones(md.numberofelements,1); +%Now, build the connectivity tables for this mesh. +md.nodeconnectivity=NodeConnectivity(md.elements,md.numberofgrids); +md.elementconnectivity=ElementConnectivity(md.elements,md.nodeconnectivity); + %outline names md.domainoutline=readfile(domainname); Index: /issm/trunk/src/m/classes/public/setelementstype.m =================================================================== --- /issm/trunk/src/m/classes/public/setelementstype.m (revision 1110) +++ /issm/trunk/src/m/classes/public/setelementstype.m (revision 1111) @@ -127,5 +127,7 @@ %icefront for stokes -md.segmentonneumann_diag_stokes=md.segmentonneumann_diag(find(md.elements_type(md.segmentonneumann_diag(:,end),2)),:); +if ~isnan(md.segmentonneumann_diag), + md.segmentonneumann_diag_stokes=md.segmentonneumann_diag(find(md.elements_type(md.segmentonneumann_diag(:,end),2)),:); +end %figure out solution types