function result_on_const_yz (i_p, y_p, z_p, out_name) % Copyright 2000, J.E. Akin. All rights reserved. % ------------------------------------------------------ % 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 BOG = -8 ; % 1st of graph tag to sepreate file saves EOG = -9 ; % end of graph tag to sepreate file saves format short g % fprintf ('Begin component value carpet plots: \n') if ( nargin == 0 ) i_p = 0 ; y_p = 0. ; z_p = 0.; out_name = ['null'] ; elseif ( nargin == 1 ) y_p = 0. ; z_p = 0.; out_name = ['null'] ; elseif ( nargin == 2 ) z_p = 0.; out_name = ['null'] ; elseif ( nargin == 3 ) out_name = ['null'] ; 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 load node_results.tmp nr = size (node_results, 1); if ( nr == 0 ) error ('Error missing file node_results.tmp') end % if error max_p = size (node_results, 2) ; % number of columns fprintf ('Read %g nodal result 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 R (np) = 0. ; % pre-allocate array R % set constants % msh_bc_xyz has: pre_p items then: x, y y = msh_bc_xyz (1:np, (pre_p+2)) ; % extract y column ymax = max (y) ; ymin = min (y) ; diff = ymax - ymin ; tol = abs ( diff / (10 * np) ) ; z = msh_bc_xyz (1:np, (pre_p+3)) ; % extract z column % create a file for combining graphs (3 items per row) %b fid = fopen('graphs_on_y.dat','a+') ; if ( out_name(1:4) ~= ['null'] ) fid = fopen(out_name,'a+') ; else fid = 0 ; end % if file given ymdhms = fix(clock) ; fprintf (fid,'%g %g %g \n', BOG, y_p, z_p ); % -8 y_p z_p fprintf (fid,'%g %g %g \n', ymdhms (1:3) ); % yr mon day fprintf (fid,'%g %g %g \n', ymdhms (4:6) ); % hr min sec given = 0 ; for j = 1:np ; % check all nodes if ( abs (y(j) - y_p) <= tol ) % found a point if ( abs (z(j) - z_p) <= tol ) % found a point given = given + 1 ; x (given) = msh_bc_xyz (j, (pre_p+1)) ; % extract y if ( i_p >= 1 ) %b if ( i_p >= 1 & i_p <=3 ) R(given) = node_results(j, i_p) ; %b elseif ( i_p == 4 ) % Von Mises Stress %b root_2 = sqrt (2.) ; %b temp = (node_results (j, 1) ... %b - node_results (j, 2))^2 ... %b + node_results (j, 1)^2 ... %b + node_results (j, 2)^2 ... %b + node_results (j, 3)^2 * 6 ; %b R (given) = sqrt ( temp ) / root_2 ; else % i_p = 0, get root mean sq R (given) = sqrt(sum(node_results(j,1:max_p).^2)); end % if get RMS value fprintf (fid,'%g %g %g \n', given, x(given), R(given) ); end % if tol on z end % if tol on y end % for j % Cite max, min values if ( given > 0 ) [V_X, L_X] = max (R(1:given)) ; [V_N, L_N] = min (R(1:given)) ; 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) if ( V_X == V_N ) % constant V_X = V_X + 2; end % if else error 'No nodes found for given position' end % if given % Initialize plots xmax = max (x(1:given)) ; xmin = min (x(1:given)) ; Rmax = max (R(1:given)) ; Rmin = min (R(1:given)) ; if ( Rmin == Rmax ) Rmax = Rmin + (xmax - xmin)/50 ; Rmin = Rmin - (xmax - xmin)/50 ; end % if % Rmin=-0.10 % Rmax=5.0 %b add delimiter to graph file for possible multiple graphs fprintf (fid,'%g %g %g \n', EOG, Rmin, Rmax ); clf % clear graphics axis ([xmin, xmax, Rmin, Rmax]) % set axes %b axis ('square') hold on % hold image for plots xlabel (['X for ', int2str(given),' Points at Y = ', ... num2str(y_p), ', Z = ', num2str(z_p), ... ' (', int2str(nt),' Elements, ', ... int2str(np),' Nodes)']) if ( i_p > 0 ) ylabel (['Component\_',int2str(i_p),' (max = ', ... num2str(V_X), ', min = ', num2str(V_N), ')']) title(['Nodal FEA Solution Component\_',int2str(i_p), ... ', (', int2str(nod_per_el), ' nodes per element)']) else % i_p = 0, get root mean sq ylabel (['Nodal RMS Solution Value (max = ', ... num2str(V_X), ', min = ', num2str(V_N), ')']) title(['Nodal RMS Flux Value']) end % if get RMS value plot(x(1:given), R(1:given),'k*') grid [y(1:given),i_sort] = sort( x(1:given) ) ; plot(x(i_sort), R(i_sort), 'k-'); %% end % if show labels % -depsc -tiff % for an eps version %b print ('-dpsc', ['result_', int2str(i_p), '_on_yz']) hold off %b v_text = ['Created result_', int2str(i_p), '_on_yz.ps'] ; %b fprintf (1,'%s', v_text) ; fprintf (1, ' \n' ) % end of result_on_const_yz