jk_calc.m

ABAQUS is a general purpose finite element analysis program which is widely used to analyses mechani

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [K_I,J_el,J_pl,Jpnd]=JK_calc(reac_set,J,num,rx,a);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% nummers  (num)						%
%     1    : ccp						%
%     2    : sec-tension					%
%     3    : sec-bending					%
%     4    : dec tension					%
%								%
% output form this file is 					%
%            J_ele   : J-elastic, calculated from K_I  		%
%            J_pl    : J-plastic, calculated as J-je		%
%	     Jpnd    : J-plastic, non dimemsionalised by	%
%		       the ligament c and yield stress s)	%
%            K_I     : K calculated from reaction forces	%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%  	Initiation of Constants Value				%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
	s0=2e8; e=2e11; po=0.3; W=1000;
	a=round(a);
	alpa=a/W;
	al=alpa;
	status=0;
	c=W-a;
%%%%%%%%%%%%%%%%%%%%%% calculation of K for J elastic %%%%%%%%%%%
%  The Stress Intensity Factor is calculated 			%
%  as described in Stress Intensity Factor Handbook, VOL 1, 	%
%  The Society of Material Science, Japan, Pergamon Press 1987	%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%								%
%%%%%%%%%%%%%% Centre Cracked Panels %%%%%%%%%%%%%%%%%%%%%%%%%%%%

	if num==1 
	   fm=[1 1.006 1.025 1.058 1.109 1.187 1.303 1.488 1.816 ...
               2.3118]; % *) 
           al=[0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9]; % *)
	   h=length(al);
             for i=1:h
                 if al(i)==alpa
                    f=fm(i);
                  end;
             end
           K_I=abs(sum(reac_set)).*(pi*a)^0.5*f/W;

%%%%%%%%%%%%%% Single Edge Tension %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

	elseif num==2
		K_I=abs(sum(reac_set)).* ...
		((pi*a)^0.5*(2*tan(pi*alpa/2)/(pi*alpa))^0.5* ...
		(0.752+2.02*alpa+0.37*(1-sin(pi*alpa/2))^3)/cos(pi*alpa/2))/W;


%%%%%%%%%%%%%% Single Edge Bend-bars %%%%%%%%%%%%%%%%%%%%%%%%%%%%
	
	elseif num==3;
	       [m,n]=size(reac_set);
		   for i=1:n;
                       M(i)=0.0;
                       for j=1:m;
                 	   M(i)= M(i)+(rx(j)-c/2)*reac_set(j,i);
                       end;
                       M=abs(M);
                   end 
	       K_I=6*M.*((pi*a)^0.5*(1.122-1.4*al+7.33*al^2 ...
		   -13.08*al^3+14*al^4)/W^2);
			
			
%%%%%%%%%%%%%% Double Edge Crack  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
			

       elseif num==4
	      	K_I=abs(sum(reac_set)).*((pi*a)^0.5* ...
	             (1.122-0.561*al-0.205*al^2+0.471*al^3 ...
	             -0.190*al^4)/(W*(1-al)^0.5));
       else 	             				                                									     	 		 	     	
                'Warning the crack problem number is not recognised'
 	         status=-1
       end;
       
%%%%%%%%%%%%%% subtraction of J elastic from J total %%%%%%%%%%%%
   if status==0
   	J_el=(K_I.^2)*(1-po^2)/e;
	J_pl=(J-J_el);
	Jpnd=J_pl/(c*s0);
   end;
end  			%( end of function) %