function color_stresses_step (i_p) % Copyright 2000, J.E. Akin. All rights reserved. % ------------------------------------------------------ % Matlab carpet plot of i_p-th component value, % at mesh node locations % If i_p = 0, show RMS value % ------------------------------------------------------ % c_x = x coordinates of nod_per_el line polygon % c_y = y coordinates of nod_per_el line polygon % msh_typ_nodes = connectivity list for elements, nt x nod_per_el % loop = corners for nod_per_el line polygon % nod_per_el = Nodes per element % np = Number of Points % nt = Number of elements % pre_e = Element items before connectivity list pre_e = 0 ; % pre_p = Nodal items before coordinates pre_p = 1; % msh_bc_xyz = Nodal coordinates (with preceeding data) % t_x = x coordinates of nod_per_el corners % t_y = y coordinates of nod_per_el corners format short g fprintf ('Begin component value carpet plots: \n') if ( nargin == 0 ) i_p = 0 ; end % if no arguments % Read coordinate file and connectivity file % integer bc code, real xy pairs for np points (pre_p = 1) load msh_bc_xyz.tmp ; % Set control data: number of points np = size (msh_bc_xyz,1) ; % number of nodal points fprintf ('Read %g mesh coordinate pairs \n', np) ns = size (msh_bc_xyz,2) - pre_p ; % space dimension if ( ns < 2 ) error ('This is not a 2D mesh') end % if not 2D data % Set control data: number elements load msh_typ_nodes.tmp ; % nod_per_el nodes per element nt = size (msh_typ_nodes,1) ; % number of elements in mesh nod_per_el = size (msh_typ_nodes,2) - pre_e -1 ; % nodes per elem fprintf ('Read %g elements connections \n', nt) load scp_node_ave_fluxes.tmp nr = size (scp_node_ave_fluxes, 1); if ( nr == 0 ) error ('Error missing file scp_node_ave_fluxes.tmp') end % if error max_p = size (scp_node_ave_fluxes, 2) ; % number of columns fprintf ('Read %g nodal stress values \n', nr) fprintf (' with %g components each \n', max_p) if ( i_p > (max_p+1) ) error ('i_p > available data') end % if error x (np) = 0. ; % pre-allocate array x y (np) = 0. ; % pre-allocate array y z (np) = 0. ; % pre-allocate array z t_nodes (nod_per_el) = 0 ; % Optional pre-allocation t_x (nod_per_el) = 0 ; % Optional pre-allocation t_y (nod_per_el) = 0 ; % Optional pre-allocation c_x (nod_per_el + 1) = 0 ; % Optional pre-allocation c_y (nod_per_el + 1) = 0 ; % Optional pre-allocation loop (nod_per_el + 1) = 0 ; % Optional pre-allocation % set constants [loop] = get_El_Loop (nod_per_el) ; % msh_bc_xyz has: pre_p items then: x, y x = msh_bc_xyz (1:np, (pre_p+1)) ; % extract x column y = msh_bc_xyz (1:np, (pre_p+2)) ; % extract y column if ( i_p >= 1 & i_p <=3 ) z = scp_node_ave_fluxes(:, i_p) ; elseif ( i_p == 4 ) % Von Mises Stress root_2 = sqrt (2.) ; for k = 1:np temp = (scp_node_ave_fluxes (k, 1) - scp_node_ave_fluxes (k, 2))^2 ... + scp_node_ave_fluxes (k, 1)^2 ... + scp_node_ave_fluxes (k, 2)^2 ... + scp_node_ave_fluxes (k, 3)^2 * 6 ; z (k) = sqrt ( temp ) / root_2 ; end % for k else % i_p = 0, get root mean sq for k = 1:np z (k) = sqrt ( sum (scp_node_ave_fluxes (k, 1:max_p).^2)) ; end % for k end % if get RMS value % Cite max, min values [V_X, L_X] = max (z) ; [V_N, L_N] = min (z) ; fprintf ('Max value is %g at node %g \n', V_X, L_X) fprintf ('Min value is %g at node %g \n', V_N, L_N) caxis([V_N V_X]) ; % set full color ranges % Initialize plots xmax = max (x) ; xmin = min (x) ; ymax = max (y) ; ymin = min (y) ; % zmax = max (z) ; zmin = min (z) ; clf % clear graphics % axis ([xmin, xmax, ymin, ymax, zmin, zmax]) % set axes axis ([xmin, xmax, ymin, ymax]) % set axes axis ('square') hold on % hold image for plots xlabel (['X at ', int2str(np),' Nodes']) ylabel (['Y on ', int2str(nt),' Elements (with ', ... int2str(nod_per_el), ' nodes) ']) if ( i_p == 1 ) %b zlabel (['Sigma\_X (max = ', ... %b num2str(V_X), ', min = ', num2str(V_N), ')']) title(['Nodal FEA SCP X-Normal Stress (max = ', ... num2str(V_X), ', min = ', num2str(V_N), ')']) elseif ( i_p == 2 ) %b zlabel (['Sigma\_Y (max = ', ... %b num2str(V_X), ', min = ', num2str(V_N), ')']) title(['Nodal FEA SCP Y-Normal Stress (max = ', ... num2str(V_X), ', min = ', num2str(V_N), ')']) elseif ( i_p == 3 ) %b zlabel (['Sigma\_XY (max = ', ... %b num2str(V_X), ', min = ', num2str(V_N), ')']) title(['Nodal FEA SCP XY-Shear Stress (max = ', ... num2str(V_X), ', min = ', num2str(V_N), ')']) elseif ( i_p == 4 ) %b zlabel (['Von Mises Criterion (max = ', ... %b num2str(V_X), ', min = ', num2str(V_N), ')']) title(['Nodal FEA SCP Von Mises Criterion (max = ', ... num2str(V_X), ', min = ', num2str(V_N), ')']) else % i_p = 0, get root mean sq %b zlabel (['Nodal RMS Stress Value (max = ', ... %b num2str(V_X), ', min = ', num2str(V_N), ')']) title(['Nodal RMS Stress Value (max = ', ... num2str(V_X), ', min = ', num2str(V_N), ')']) end % if get RMS value % Loop over all elements for it = 1:nt ; % Extract corner connectivity t_nodes = msh_typ_nodes (it, (pre_e+2):(nod_per_el+pre_e+1)); % Extract corner coordinates t_x = x (t_nodes) ; % x at those nodes, only t_y = y (t_nodes) ; % y at those nodes, only t_z = z (t_nodes) ; % z at those nodes, only % Plot this polygon c_x = t_x (loop) ; % x for nod_per_el line polygon c_y = t_y (loop) ; % y for nod_per_el line polygon c_z = t_z (loop) ; color = sum (c_z' )/nod_per_el ; fill (c_x, c_y, color) % plot nod_per_el lines end % for over all elements fill (c_x, c_y, c_z), grid % add grid to last one colorbar % label max min points text ('Color', 'yellow') v_text = sprintf ('------min') ; text (x(L_N), y(L_N), V_N, v_text) ; v_text = sprintf ('------max') ; text (x(L_X), y(L_X), V_X, [v_text]) ; ; % end % if show labels % -depsc -tiff % for an eps version %bprint ('-dpsc', ['color_stresses_', int2str(i_p), '_step']) hold off %bv_text = ['Created color_stresses_', int2str(i_p), '_step.ps'] ; %bfprintf (1,'%s', v_text) ; fprintf (1, ' \n' ) % end of color_stresses_step