[cig-commits] r5611 - in long/3D/Gale/trunk: . tools
walter at geodynamics.org
walter at geodynamics.org
Thu Dec 21 00:26:31 PST 2006
Author: walter
Date: 2006-12-21 00:26:31 -0800 (Thu, 21 Dec 2006)
New Revision: 5611
Added:
long/3D/Gale/trunk/tools/PlotGaleOutput.m
long/3D/Gale/trunk/tools/ReadGaleOutput.m
long/3D/Gale/trunk/tools/ReshapeGaleData.m
Modified:
long/3D/Gale/trunk/
Log:
r1257 at earth: boo | 2006-12-21 00:24:14 -0800
Add matlab vis files
Property changes on: long/3D/Gale/trunk
___________________________________________________________________
Name: svk:merge
- 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:1256
+ 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:1257
Added: long/3D/Gale/trunk/tools/PlotGaleOutput.m
===================================================================
--- long/3D/Gale/trunk/tools/PlotGaleOutput.m 2006-12-21 08:26:27 UTC (rev 5610)
+++ long/3D/Gale/trunk/tools/PlotGaleOutput.m 2006-12-21 08:26:31 UTC (rev 5611)
@@ -0,0 +1,155 @@
+%PlotGaleOutput
+%
+% Reads in Gale data and creates plots
+
+
+clear, close all
+
+fname = 'extension'; % directory-name
+step = 1 ; % timestep
+
+
+% Read data
+ReadGaleOutput;
+
+% Color definitions etc.
+MeshColor = 'k';
+VelocityArrowColor = 'k';
+
+v=version;
+Octave=v(1)<'5';
+
+if (Octave)
+ toggle_octplot
+end
+
+% Plot data
+if (dim==2)
+
+
+ if (Octave)
+ figure, clf
+ else
+ figure(111), clf
+ quiver(Velocity.x,Velocity.y,Velocity.Vx,Velocity.Vy,VelocityArrowColor); % Velocity arrows
+ end
+ hold on
+ %pcolor(Velocity.x,Velocity.y,Velocity.Vx); % Background color: Vx velocity
+ pcolor(StrainRateInvariant.x,StrainRateInvariant.y,StrainRateInvariant.data); % Background color: Vx velocity
+
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ if (Octave==0)
+ axis('equal')
+ axis('tight')
+ end
+ title(['\epsilon_{II}; time=',num2str(time),' dt=',num2str(dt),' step=',num2str(step)])
+
+ if(Octave)
+ figure
+ else
+ figure(112)
+ end
+
+ subplot(331)
+ pcolor(StrainRateInvariant.x,StrainRateInvariant.y,StrainRateInvariant.data);
+ hold on
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ if(Octave==0)
+ colorbar
+ end
+ title('StrainRateInvariant');
+
+ subplot(332)
+ pcolor(StrainRate.x,StrainRate.y,StrainRate.data);
+ hold on
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ if(Octave==0)
+ colorbar
+ end
+ title('StrainRate');
+
+ subplot(333)
+ pcolor(Vorticity.x,Vorticity.y,Vorticity.xy);
+ hold on
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ if(Octave==0)
+ colorbar
+ end
+ title('Vorticity_{xy}');
+
+ subplot(334)
+ pcolor(Pressure.x,Pressure.y,Pressure.data);
+ hold on
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ if(Octave==0)
+ colorbar
+ end
+ title('Pressure');
+
+ subplot(335)
+ pcolor(StrainRateXX.x,StrainRateXX.y,StrainRateXX.data);
+ hold on
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ if(Octave==0)
+ colorbar
+ end
+ title('StrainRateXX');
+
+ subplot(336)
+ pcolor(StrainRateYY.x,StrainRateYY.y,StrainRateYY.data);
+ hold on
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ if(Octave==0)
+ colorbar
+ end
+ title('StrainRateYY');
+
+ subplot(337)
+ pcolor(VelocityGradientsField.x,VelocityGradientsField.y,VelocityGradientsField.data);
+ hold on
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ if(Octave==0)
+ colorbar
+ end
+ title('VelocityGradientsField');
+
+ subplot(338)
+ pcolor(VelocityMagnitudeField.x,VelocityMagnitudeField.y,VelocityMagnitudeField.data);
+ hold on
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ if(Octave==0)
+ colorbar
+ end
+ title('VelocityMagnitudeField');
+
+ if(Octave)
+ figure
+ else
+ subplot(339)
+ end
+ plot(Yielding.x,Yielding.y,'.');
+ hold on
+ plot(mesh_x,mesh_y,MeshColor,mesh_x.',mesh_y.',MeshColor); % FEM mesh
+ title('Particles');
+ if(Octave==0)
+ axis('tight')
+ end
+
+elseif (dim==3)
+
+ figure(1), clf
+ slice(StrainRateInvariant.x,StrainRateInvariant.z,StrainRateInvariant.y,StrainRateInvariant.data,[0 1],[1],[.1 ]), hold on
+ quiver3(Velocity.x(:),Velocity.z(:),Velocity.y(:),Velocity.Vx(:),Velocity.Vz(:),Velocity.Vy(:),VelocityArrowColor); % Velocity arrows
+ axis('equal')
+ if(Octave==0)
+ box('on')
+ end
+ title('\epsilon_{II} and velocity (arrows)');
+ if(Octave==0)
+ colorbar
+ end
+
+
+
+
+end
Added: long/3D/Gale/trunk/tools/ReadGaleOutput.m
===================================================================
--- long/3D/Gale/trunk/tools/ReadGaleOutput.m 2006-12-21 08:26:27 UTC (rev 5610)
+++ long/3D/Gale/trunk/tools/ReadGaleOutput.m 2006-12-21 08:26:31 UTC (rev 5611)
@@ -0,0 +1,192 @@
+%ReadGaleOutput
+%
+% Reads a timestep with Gale output data
+
+
+% To be done before release:
+% - check 3D case, and generation of 3D matrixes
+
+
+%fname = 'extension';
+%step = 2;
+
+addpath(pwd); %adds current directory to path
+
+
+% Change to output directory
+dirname = ['../output.',fname];
+cd(dirname);
+
+% Read the input xml-file, to get the number of nodes in each direction----
+%
+fid = fopen('input.xml');
+while ~feof(fid)
+ line = fgetl(fid);
+
+ if(ischar(line)==0)
+ break
+ end
+ in = strfind(line,'"elementResI"');
+ if ~isempty(in);
+ indstart = strfind(line,'>');
+ indend = strfind(line,'<');
+ nelI = str2num( line(indstart(1)+1:indend(end)-1 ) );
+ end
+
+ in = strfind(line,'"elementResJ"');
+ if ~isempty(in);
+ indstart = strfind(line,'>');
+ indend = strfind(line,'<');
+ nelJ = str2num( line(indstart(1)+1:indend(end)-1 ) );
+ end
+
+ in = strfind(line,'"elementResK"');
+ if ~isempty(in);
+ indstart = strfind(line,'>');
+ indend = strfind(line,'<');
+ nelK = str2num( line(indstart(1)+1:indend(end)-1 ) );
+ end
+
+ in = strfind(line,'"dim"');
+ if ~isempty(in);
+ indstart = strfind(line,'>');
+ indend = strfind(line,'<');
+ dim = str2num( line(indstart(1)+1:indend(end)-1 ) );
+ end
+
+end
+fclose(fid);
+
+if dim==2
+ NumberElements = [nelI nelJ ];
+ NumberNodes = [nelI+1 nelJ+1 ];
+
+elseif dim==3
+ NumberElements = [nelI nelJ nelK ];
+ NumberNodes = [nelI+1 nelJ+1 nelK+1 ];
+
+end
+
+%
+%--------------------------------------------------------------------------
+
+
+% Read all properties -----------------------------------------------------
+%
+
+% Load the data
+number=num2str(step+100000,'%5.0f'); number(1)=[];
+
+% Check that file-number exists
+if exist(['VelocityField','.',number,'.dat'])
+
+ % Pressure (constant per element)
+ data = load(['PressureField','.',number,'.dat']);
+ Pressure.x = ReshapeGaleData(dim,NumberElements,data(:,2) );
+ Pressure.y = ReshapeGaleData(dim,NumberElements,data(:,3) );
+ Pressure.data = ReshapeGaleData(dim,NumberElements,data(:,5) );
+ if dim==3
+ Pressure.z = ReshapeGaleData(dim,NumberElements,data(:,4) );
+ end
+
+ % Velocity
+ data = load(['VelocityField','.',number,'.dat']);
+ Velocity.x = ReshapeGaleData(dim,NumberNodes,data(:,2) );
+ Velocity.y = ReshapeGaleData(dim,NumberNodes,data(:,3) );
+ Velocity.z = ReshapeGaleData(dim,NumberNodes,data(:,4) );
+
+ Velocity.Vx = ReshapeGaleData(dim,NumberNodes,data(:,5) );
+ Velocity.Vy = ReshapeGaleData(dim,NumberNodes,data(:,6) );
+ if dim==3
+ Velocity.Vz = ReshapeGaleData(dim,NumberNodes,data(:,7) );
+ end
+
+ % StrainRateInvariantField
+ data = load(['StrainRateInvariantField','.',number,'.dat']);
+ StrainRateInvariant.x = ReshapeGaleData(dim,NumberNodes,data(:,2) );
+ StrainRateInvariant.y = ReshapeGaleData(dim,NumberNodes,data(:,3) );
+ StrainRateInvariant.z = ReshapeGaleData(dim,NumberNodes,data(:,4) );
+ StrainRateInvariant.data = ReshapeGaleData(dim,NumberNodes,data(:,5) );
+
+ % VorticityField
+ data = load(['VorticityField','.',number,'.dat']);
+ Vorticity.x = ReshapeGaleData(dim,NumberNodes,data(:,2) );
+ Vorticity.y = ReshapeGaleData(dim,NumberNodes,data(:,3) );
+ Vorticity.z = ReshapeGaleData(dim,NumberNodes,data(:,4) );
+ Vorticity.xx = ReshapeGaleData(dim,NumberNodes,data(:,5) );
+ Vorticity.xy = ReshapeGaleData(dim,NumberNodes,data(:,6) );
+ Vorticity.yx = ReshapeGaleData(dim,NumberNodes,data(:,7) );
+ Vorticity.yy = ReshapeGaleData(dim,NumberNodes,data(:,8) );
+
+
+ % StrainRateField
+ data = load(['StrainRateField','.',number,'.dat']);
+ StrainRate.x = ReshapeGaleData(dim,NumberNodes,data(:,2) );
+ StrainRate.y = ReshapeGaleData(dim,NumberNodes,data(:,3) );
+ StrainRate.z = ReshapeGaleData(dim,NumberNodes,data(:,4) );
+ StrainRate.data = ReshapeGaleData(dim,NumberNodes,data(:,5) );
+
+ % StrainRateXXField
+ data = load(['StrainRateXXField','.',number,'.dat']);
+ StrainRateXX.x = ReshapeGaleData(dim,NumberNodes,data(:,2) );
+ StrainRateXX.y = ReshapeGaleData(dim,NumberNodes,data(:,3) );
+ StrainRateXX.z = ReshapeGaleData(dim,NumberNodes,data(:,4) );
+ StrainRateXX.data = ReshapeGaleData(dim,NumberNodes,data(:,5) );
+
+ % StrainRateYYField
+ data = load(['StrainRateYYField','.',number,'.dat']);
+ StrainRateYY.x = ReshapeGaleData(dim,NumberNodes,data(:,2) );
+ StrainRateYY.y = ReshapeGaleData(dim,NumberNodes,data(:,3) );
+ StrainRateYY.z = ReshapeGaleData(dim,NumberNodes,data(:,4) );
+ StrainRateYY.data = ReshapeGaleData(dim,NumberNodes,data(:,5) );
+
+ % StrainRateYYField
+ data = load(['StrainRateYYField','.',number,'.dat']);
+ StrainRateYY.x = ReshapeGaleData(dim,NumberNodes,data(:,2) );
+ StrainRateYY.y = ReshapeGaleData(dim,NumberNodes,data(:,3) );
+ StrainRateYY.z = ReshapeGaleData(dim,NumberNodes,data(:,4) );
+ StrainRateYY.data = ReshapeGaleData(dim,NumberNodes,data(:,5) );
+
+ % VelocityGradientsField
+ data = load(['VelocityGradientsField','.',number,'.dat']);
+ VelocityGradientsField.x = ReshapeGaleData(dim,NumberNodes,data(:,2) );
+ VelocityGradientsField.y = ReshapeGaleData(dim,NumberNodes,data(:,3) );
+ VelocityGradientsField.z = ReshapeGaleData(dim,NumberNodes,data(:,4) );
+ VelocityGradientsField.data = ReshapeGaleData(dim,NumberNodes,data(:,5) );
+
+ % VelocityMagnitudeField
+ data = load(['VelocityMagnitudeField','.',number,'.dat']);
+ VelocityMagnitudeField.x = ReshapeGaleData(dim,NumberNodes,data(:,2) );
+ VelocityMagnitudeField.y = ReshapeGaleData(dim,NumberNodes,data(:,3) );
+ VelocityMagnitudeField.z = ReshapeGaleData(dim,NumberNodes,data(:,4) );
+ VelocityMagnitudeField.data = ReshapeGaleData(dim,NumberNodes,data(:,5) );
+
+
+ % Yielding
+ data = load(['yielding','.',number,'.txt']);
+ Yielding.x = data(:,2);
+ Yielding.y = data(:,3);
+ Yielding.z = data(:,4);
+ Yielding.data = data(:,5);
+
+ % Time_info
+ data = load(['timeInfo','.',number,'.dat']);
+ time = data(1);
+ dt = data(2);
+
+ % Mesh info
+ mesh_x = Velocity.x;
+ mesh_y = Velocity.y;
+
+
+else
+ cd ..
+
+ error(['FileNumber ',number,' does not exist!'])
+
+end
+%
+%--------------------------------------------------------------------------
+
+% Change back to previous directory
+cd ..
Added: long/3D/Gale/trunk/tools/ReshapeGaleData.m
===================================================================
--- long/3D/Gale/trunk/tools/ReshapeGaleData.m 2006-12-21 08:26:27 UTC (rev 5610)
+++ long/3D/Gale/trunk/tools/ReshapeGaleData.m 2006-12-21 08:26:31 UTC (rev 5611)
@@ -0,0 +1,12 @@
+function [data] = ReshapeGaleData(dim, NumberProps, data_vector);
+% Reshapes Gale data into 2D or 3D arrays
+%
+
+
+if (dim==2)
+ data = reshape(data_vector, NumberProps);
+elseif (dim==3)
+
+ data = reshape(data_vector, NumberProps);
+ data = shiftdim(data,2);
+end
More information about the cig-commits
mailing list