Changeset 3873
- Timestamp:
- 05/21/10 11:04:58 (15 years ago)
- Location:
- issm/trunk/src/m/solutions/jpl
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk/src/m/solutions/jpl/ModelUpdateInputsFromVector.m
r3853 r3873 1 function ModelUpdateInputsFromVector(models, vector, enum, typeenum);1 function models=ModelUpdateInputsFromVector(models, vector, enum, typeenum); 2 2 %MODELUPDATEINPUTSFROMVECTOR - Update inputs using a vector 3 3 % … … 26 26 if isstruct(model), %there is an analysis_type model 27 27 [model.elements,model.nodes,model.vertices,model.loads,model.materials,model.parameters] = UpdateInputsFromVector(model.elements,model.nodes,model.vertices,model.loads,model.materials,model.parameters,vector,enum, typeenum); 28 models.(field)=model; 28 29 end 29 30 -
issm/trunk/src/m/solutions/jpl/diagnostic_core.m
r3853 r3873 6 6 % 7 7 8 %recover fem models9 m_dh=model.dh;10 m_dv=model.dv;11 m_ds=model.ds;12 m_dhu=model.dhu;13 m_sl=model.sl;14 15 8 %recover parameters common to all solutions 16 verbose=m _dhu.parameters.Verbose;17 dim=m _dh.parameters.Dim;18 ishutter=m _dhu.parameters.IsHutter;19 ismacayealpattyn=m _dh.parameters.IsMacAyealPattyn;20 isstokes=m _ds.parameters.IsStokes;21 numrifts=m _dhu.parameters.NumRifts;22 qmu_analysis=m _dh.parameters.QmuAnalysis;9 verbose=model.dhu.parameters.Verbose; 10 dim=model.dh.parameters.Dim; 11 ishutter=model.dhu.parameters.IsHutter; 12 ismacayealpattyn=model.dh.parameters.IsMacAyealPattyn; 13 isstokes=model.ds.parameters.IsStokes; 14 numrifts=model.dhu.parameters.NumRifts; 15 qmu_analysis=model.dh.parameters.QmuAnalysis; 23 16 24 17 25 18 %for qmu analysis, be sure the velocity input we are starting from is the one in the parameters: 26 19 if qmu_analysis, 27 ModelUpdateInputsFromVector(model,m_dh.vx,VxEnum,VertexEnum);28 ModelUpdateInputsFromVector(model,m_dh.vy,VyEnum,VertexEnum);29 ModelUpdateInputsFromVector(model,m_dh.vz,VzEnum,VertexEnum);20 model=ModelUpdateInputsFromVector(model,m_dh.vx,VxEnum,VertexEnum); 21 model=ModelUpdateInputsFromVector(model,m_dh.vy,VyEnum,VertexEnum); 22 model=ModelUpdateInputsFromVector(model,m_dh.vz,VzEnum,VertexEnum); 30 23 end 31 24 32 25 %Compute slopes: 33 [surfaceslopex,surfaceslopey]=slope_core(m _sl,SurfaceAnalysisEnum);34 [bedslopex,bedslopey]=slope_core(m _sl,BedAnalysisEnum);26 [surfaceslopex,surfaceslopey]=slope_core(model.sl,SurfaceAnalysisEnum); 27 [bedslopex,bedslopey]=slope_core(model.sl,BedAnalysisEnum); 35 28 36 29 %Update: 37 ModelUpdateInputsFromVector(model,surfaceslopex,SurfaceSlopexEnum,VertexEnum);38 ModelUpdateInputsFromVector(model,surfaceslopey,SurfaceSlopeyEnum,VertexEnum);39 ModelUpdateInputsFromVector(model,bedslopex,BedSlopexEnum,VertexEnum);40 ModelUpdateInputsFromVector(model,bedslopey,BedSlopeyEnum,VertexEnum);30 model=ModelUpdateInputsFromVector(model,surfaceslopex,SurfaceSlopexEnum,VertexEnum); 31 model=ModelUpdateInputsFromVector(model,surfaceslopey,SurfaceSlopeyEnum,VertexEnum); 32 model=ModelUpdateInputsFromVector(model,bedslopex,BedSlopexEnum,VertexEnum); 33 model=ModelUpdateInputsFromVector(model,bedslopey,BedSlopeyEnum,VertexEnum); 41 34 42 35 if ishutter, 43 36 44 37 displaystring(verbose,'\n%s',['computing hutter velocities...']); 45 u_g=diagnostic_core_linear(m _dhu,DiagnosticAnalysisEnum(),HutterAnalysisEnum());38 u_g=diagnostic_core_linear(model.dhu,DiagnosticAnalysisEnum(),HutterAnalysisEnum()); 46 39 47 40 displaystring(verbose,'\n%s',['computing pressure according to MacAyeal...']); 48 p_g=ComputePressure(m _dhu.elements,m_dhu.nodes,m_dhu.vertices,m_dhu.loads,m_dhu.materials,m_dhu.parameters,DiagnosticAnalysisEnum(),HutterAnalysisEnum());41 p_g=ComputePressure(model.dhu.elements,model.dhu.nodes,model.dhu.vertices,model.dhu.loads,model.dhu.materials,model.dhu.parameters,DiagnosticAnalysisEnum(),HutterAnalysisEnum()); 49 42 50 43 displaystring(verbose,'\n%s',['update boundary conditions for macyeal pattyn using hutter results...']); 51 44 if ismacayealpattyn, 52 m _dh.y_g=u_g;53 [m _dh.ys m_dh.ys0]=Reducevectorgtos(m_dh.y_g,m_dh.nodesets);45 model.dh.y_g=u_g; 46 [model.dh.ys model.dh.ys0]=Reducevectorgtos(model.dh.y_g,model.dh.nodesets); 54 47 end 55 48 … … 59 52 60 53 displaystring(verbose,'\n%s',['computing horizontal velocities...']); 61 [u_g m _dh.loads]=diagnostic_core_nonlinear(m_dh,m_dh.loads,DiagnosticAnalysisEnum(),HorizAnalysisEnum());54 [u_g model.dh.loads]=diagnostic_core_nonlinear(model.dh,model.dh.loads,DiagnosticAnalysisEnum(),HorizAnalysisEnum()); 62 55 63 56 if dim==2, 64 57 displaystring(verbose,'\n%s',['computing pressure according to MacAyeal...']); 65 p_g=ComputePressure(m _dh.elements,m_dh.nodes,m_dh.vertices,m_dh.loads,m_dh.materials,m_dh.parameters,DiagnosticAnalysisEnum(),HorizAnalysisEnum());58 p_g=ComputePressure(model.dh.elements,model.dh.nodes,model.dh.vertices,model.dh.loads,model.dh.materials,model.dh.parameters,DiagnosticAnalysisEnum(),HorizAnalysisEnum()); 66 59 end 67 60 end … … 70 63 71 64 displaystring(verbose,'\n%s',['extruding horizontal velocities...']); 72 u_g_horiz=FieldExtrude(m _dh.elements,m_dh.nodes,m_dh.vertices,m_dh.loads,m_dh.materials,m_dh.parameters,u_g,'velocity',1);65 u_g_horiz=FieldExtrude(model.dh.elements,model.dh.nodes,model.dh.vertices,model.dh.loads,model.dh.materials,model.dh.parameters,u_g,'velocity',1); 73 66 74 [vx,vy]=SplitSolutionVector(u_g_horiz,m _dh.parameters.NumberOfNodes,m_dh.parameters.NumberOfDofsPerNode);75 ModelUpdateInputsFromVector(model,vx,VxEnum,VertexEnum);76 ModelUpdateInputsFromVector(model,vy,VyEnum,VertexEnum);67 [vx,vy]=SplitSolutionVector(u_g_horiz,model.dh.parameters.NumberOfNodes,model.dh.parameters.NumberOfDofsPerNode); 68 model=ModelUpdateInputsFromVector(model,vx,VxEnum,VertexEnum); 69 model=ModelUpdateInputsFromVector(model,vy,VyEnum,VertexEnum); 77 70 78 71 displaystring(verbose,'\n%s',['computing vertical velocities...']); 79 u_g_vert=diagnostic_core_linear(m _dv,DiagnosticAnalysisEnum(),VertAnalysisEnum());80 ModelUpdateInputsFromVector(model,u_g_vert,VzEnum,VertexEnum);72 u_g_vert=diagnostic_core_linear(model.dv,DiagnosticAnalysisEnum(),VertAnalysisEnum()); 73 model=ModelUpdateInputsFromVector(model,u_g_vert,VzEnum,VertexEnum); 81 74 82 75 displaystring(verbose,'\n%s',['computing pressure according to Pattyn...']); 83 p_g=ComputePressure(m_dh.elements,m_dh.nodes,m_dh.vertices,m_dh.loads,m_dh.materials,m_dh.parameters,DiagnosticAnalysisEnum(),HorizAnalysisEnum()); 84 ModelUpdateInputsFromVector(model,p_g,PressureEnum,VertexEnum); 76 p_g=ComputePressure(model.dh.elements,model.dh.nodes,model.dh.vertices,model.dh.loads,model.dh.materials,model.dh.parameters,DiagnosticAnalysisEnum(),HorizAnalysisEnum()); 77 model=ModelUpdateInputsFromVector(model,p_g,PressureEnum,VertexEnum); 78 u_g=zeros(model.dh.parameters.NumberOfNodes*3,1); %%%%%%%%%%%%%%%%%%%%%%% NEED TO BE CLEANED 79 u_g(1:3:end)=vx; 80 u_g(2:3:end)=vy; 81 u_g(3:3:end)=u_g_vert; 85 82 86 83 if isstokes, … … 90 87 91 88 displaystring(verbose,'\n%s',['update boundary conditions for stokes using velocities previously computed...']); 92 m _ds.y_g=zeros(m_ds.nodesets.gsize,1);93 m _ds.y_g(dofsetgen([1,2],4,m_ds.nodesets.gsize))=u_g;94 m _ds.y_g(dofsetgen([3],4,m_ds.nodesets.gsize))=u_g_vert;95 [m _ds.ys m_ds.ys0]=Reducevectorgtos(m_ds.y_g,m_ds.nodesets);89 model.ds.y_g=zeros(model.ds.nodesets.gsize,1); 90 model.ds.y_g(dofsetgen([1,2],4,model.ds.nodesets.gsize))=u_g; 91 model.ds.y_g(dofsetgen([3],4,model.ds.nodesets.gsize))=u_g_vert; 92 [model.ds.ys model.ds.ys0]=Reducevectorgtos(model.ds.y_g,model.ds.nodesets); 96 93 97 94 displaystring(verbose,'\n%s',['computing stokes velocities and pressure ...']); 98 u_g=diagnostic_core_nonlinear(m _ds,DiagnosticAnalysisEnum(),StokesAnalysisEnum());95 u_g=diagnostic_core_nonlinear(model.ds,DiagnosticAnalysisEnum(),StokesAnalysisEnum()); 99 96 100 97 %"decondition" pressure 101 p_g=u_g(4:4:end)*m _dh.parameters.stokesreconditioning;98 p_g=u_g(4:4:end)*model.dh.parameters.stokesreconditioning; 102 99 end 103 100 end … … 109 106 110 107 if numrifts, 111 results.riftproperties=OutputRifts(m _dh.loads,m_dh.parameters);108 results.riftproperties=OutputRifts(model.dh.loads,model.dh.parameters); 112 109 end -
issm/trunk/src/m/solutions/jpl/diagnostic_core_linear.m
r3839 r3873 8 8 m.parameters.Kflag=1; m.parameters.Pflag=1; 9 9 10 %Update inputs in datasets11 [m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters);12 13 10 %system matrices 14 11 [K_gg, p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,analysis_type,sub_analysis_type); … … 21 18 %Reduce load from g size to f size 22 19 [p_f] = Reduceloadfromgtof( p_g, m.Gmn, K_fs, m.ys, m.nodesets); 23 20 24 21 %Solve 25 22 [u_f]=Solver(K_ff,p_f,[],m.parameters);
Note:
See TracChangeset
for help on using the changeset viewer.