Solution

cyl2car.m

function [xc,yc]=cyl2car(r,thet) %function 6: conversion of coordinates

xc=r*cos(thet);

yc=r*sin(thet);

end

function lf1=loadf1(estress,ystress) %function 4: Load Factor I

lf1=estress/ystress;

end

function lf2=loadf2(mstress,c,ftough) %function 5: Load Factor II

lf2=mstress*sqrt(c)/ftough;

end

%Declaration of parameters and other constants

clear all

clc

opengl software

c=0.1; %diameter of circular crack

P=10E3; %internal pressure

T=2E6; %torque

sigax=10E3; %applied axial stress

M1=-1E5; %bending moment along Y axis

M2=1E6; %bending moment along X axis

ystress=68200; %yield stress

ftough=7.252; %fracture toughness

%mainroutine: calculates and plots the contours for the load factors

r=a:(b-a)/1000:b; %radial coordinate array

thet=0:(2*pi)/1000:2*pi; %angular coordinate array

[sr1,sr2]=size(r);

[st1,st2]=size(thet);

xxm=zeros(sr2,st2);

yym=zeros(sr2,st2);

LF1m=zeros(sr2,st2);

LF2m=zeros(sr2,st2);

for j=1:1:sr2

for k=1:1:st2

[xx(j,k),yy(j,k)]=cyl2car(r(j),thet(k));

[sigx,sigy,sigapp]=nstress(M1,M2,xx(j,k),yy(j,k),M1+M2,sigax);

tt=tstress(T,r(j),M1+M2);

end

end

figure(1) %start of source code

axes(‘FontSize’,18)

hold on

colormap(jet)

axis square

pcolor(xx,yy,LF1)

numcount=5;

contourf(xx,yy,LF1),xlabel(‘x position (in)’),ylabel(‘y position (in)’),title(‘LF1’)

colorbar(‘southoutside’)

hold off

figure(2)

axes(‘FontSize’,18)

hold on

colormap(jet)

axis square

pcolor(xx,yy,LF2)

numcount=5;

contourf(xx,yy,LF2),xlabel(‘x position (in)’),ylabel(‘y position (in)’),title(‘LF2’)

colorbar(‘southoutside’)

hold off  %end of source code

nstress.m

function [sigxsigysigapp]=nstress(M1,M2,x,y,I,sigax) %function 2:stresses due to applied axial stress, and bending moments

sigx=(M1*x/I);

sigy=(M2*y/I);

sigapp=sigax;

end

pstress.m

function [sigrad,sigtan]=pstress(P,a,b,r) %function 1: stresses due to pressure

sigrad=-P*(b^2)*((a^2)-(r^2))/((r^2)*((a^2)-(b^2))); %Stress due to pressure in radial direction

sigtan=P*(b^2)*((a^2)+(r^2))/((r^2)*((a^2)-(b^2))); %Stress due to pressure in tangential direction

end

tstress.m

functiontt=tstress(T,r,J) %function 3: stress due to torque

tt=(T*r/J);

end