|
duci
Forum Admin
172 Posts |
 Posted - 08/07/2005 : 10:15:40
|
program polfig;
{Program pentru plotarea cu linii de egal nivel.
Funtia pentru reprezentare grafica este data de
variatia spatiala a modulu lui de elasticitate
pentru Zircaloy-4, un aliaj de Zr cu aplicatii
în energetica nucleara, caracterizat de o
"textura" accentuata }
uses crt, graph;
type vector=array [1..640] of real;
matrice=array [0..32,0..32] of real;
const { Coeficientii de textura}
{ ------- Recristalizat alfa-------- }
c000 = 1;
c200 = 2.550;
c220 = 0.493;
c400 = 0.868;
c420 = 1.070;
c440 = 0.149;
var l,i,j,m,m9,ix,iy,jx,jy :integer;
scx, scy, xp, yp, xx, yy :real;
minim,maxim,n00,nx1,ny1,n11 :real;
nivel,nivel0,t,tp,t3,temp :real;
fi,alfa,x,y,izotrop :real;
stri,stri1 :string;
a :matrice;
{ ------ Coeficientii de elasticitate ---------}
s11,s13,s12,s33,s44 :real;
{~ Modulii elastici E,G si E^-1, G^-1 ~}
e,g,em1,gm1 : real;
aa,bb,cc,ae,be,ce,ag,bg,cg : real;
{~----------- Definirea functiilor sferice de suprafata --------~}
function k20 (fi,alfa:real) :real;
Begin
k20 := 1/4*sqrt(5/pi)*(3*sqr(cos(fi))-1);
End;
function k22 (fi,alfa:real) :real;
Begin
k22 := -1/4*sqrt(15/2/pi)*sqr(sin(fi))*(cos(2*alfa));
End;
function k40 (fi,alfa:real) :real;
var cp :real;
Begin
cp := sqr(cos(fi));
k40 := 6/16/sqrt(2*pi)*(35*sqr(cp)-30*cp+3);
End;
function k42(fi,alfa:real):real;
var cp,sp:real;
Begin
cp := sqr(cos(fi));
sp := 1-cp;
k42 := -6/16*sqrt(5/2/pi)*sp*(7*cp-1)*cos(2*alfa);
End;
function k44(fi,alfa:real) :real;
var sp :real;
Begin
sp := sqr(sin(fi));
k44 := 3/16*sqrt(35/2/pi)*sqr(sp)*(cos(4*alfa));
End;
function max(a,b:real) :real;
Begin
if a>b then max:=a else max:=b;
End;
function min(a,b:real) :real;
Begin
if a<b then min:=a else min:=b;
End;
{ ---Transformarea din coordonate carteziene in
coordonatele sferice pentru proiectia stereografica -----}
procedure transformare(ax,ay:real;var ix,iy:integer);
var fi,alfa,x,y :real;
Begin
fi := ax*pi/120;{ Fi/2 in radiani}
alfa := 3*ay*pi/180;
x := 30*sin(fi)/cos (fi) ; {Raza in cercul de proiectie}
y := x*sin(alfa) ;{ Coordonatele de plotare}
x := x*cos(alfa) ;
ix := 40+round((x)*scx) ;
iy := 320-round((y)*scy) ;
End;
procedure sus;
var fi, alfa,x,y :real;
ix,iy :integer;
Begin
if n00*nx1<0 then
begin
transformare(i+n00/(n00-nx1),j,ix,iy);
lineto(ix,iy);
end;
End;
{------- Testarea si plotarea liniilor de egal nivel
pe diferitele laturi ale celulei luata in calcul----}
procedure stanga;
var fi,alfa,x,y: real;
ix,iy :integer;
Begin
if n00*ny1<0 then
begin
transformare(i,j+n00/(n00-ny1),ix,iy);
lineto(ix,iy);
end;
End;
procedure jos;
var fi, alfa,x,y :real;
ix,iy :integer;
Begin
if n11*ny1<0 then
begin
transformare(i+ny1/(ny1-n11),j+1,ix,iy);
lineto(ix,iy);
end;
End;
procedure dreapta;
var fi,alfa,x,y :real;
ix,iy :integer;
Begin
if n11*nx1<0 then
begin
transformare(i+1,j+nx1/(nx1-n11),ix,iy);
lineto(ix,iy);
end;
End;
{-------------- Proramul principal --------------}
BEGIN
t := 20;{ temperatura in grade Celsius }
tp := sqr(t);
t3 := t* tp;
s11 := 9.855+0.0104*t-0.61e-5*tp+2e-8*t3; {=E^-1 }
s13 := -2.391+9.11e-4*t+1.1e-6*tp;
s12 := -3.797+9.27e-3*t+7.51e-6*tp+2e-8*t3;
s33 := 7.925+2.26e-3*t+1.52e-6*tp;
s44 := 30.874+0.0157*t+4.39e-6*tp+1e-8*t3;
{
E:=ae*sin^4(teta)+be*cos^4(teta)+ce*sin^2(teta)*cos^2(teta)
G:=ag+bg*sin^2(teta)+bg*cos^2(teta)*sin^2(teta)
}
temp := s33*(s11+s12) -2*sqr(s13);
ag := 1/s44;
bg := 1/(s11-s12)-1/2/s44;
cg := (s33+2*s11+2*s12+4*s13)/temp+1/(s11-s12)-2/s44;
aa := ag+bg;
bb := cg-bg;
cc := -cg;
{---- Calculul modulilor elastici pentru proba texturata ----}
izotrop := (aa+bb/3+cc/5);
for i:= 0 to 31 do
begin
for j := 0 to 31 do
begin
fi := (3*i)*pi/180;
alfa := (3*j)*pi/180;
temp := izotrop+
sqrt(pi/125)*(40/35*cc+4*bb/3)*
(c200*k20(fi,alfa)+c220*k20(fi,alfa))
+8/9*cc/105*sqrt(2*pi)*
(c400*k40(fi,alfa)+c420*k42(fi,alfa)+
c440*k44(fi, alfa));
a[i,j]:=temp ;
end;
end;
{-------------------- Scalarea datelor -----------------------}
minim := a[0,0];
maxim := minim;
for i:= 0 to 30 do
begin
for j := 0 to 30 do
begin
minim := min(minim, a[i,j]);
maxim := max(maxim, a[i,j]);
end;
end;
m9 := 19;{ Nr. claselor de egal nivel }
scy := 320/30;//240/(30);{ Factorul de scala pe axa y }
scx := 320/(30); { Factorul de scala pe axa x }
nivel0 := (maxim-minim)/m9;
nivel := minim;
i := detect;j := detect;
initgraph(i,j,'d:\tp\bgi');
cleardevice;
l:=getbkcolor;
outtextxy (410,1,' Textura "Recristalizat alfa"');
outtextxy (410,21, ' Modulul de forfecare/1000 [GPa] , 20 C');
{------ Trasarea cercului pentru proiectia stereografica -----}
line (39,1,39,321);
line (39,321,360,321);
for i := 0 to 9 do
begin
fi := (i)*pi/36;
x := sin(fi) /cos (fi);
ix := 40+round(x*scx*30);
iy := 320-round(x*scy*30);
line(30,iy,39,iy);
str(i*10,stri);
outtextxy(15,iy-4,stri);
line (ix, 321, ix,331) ;
str(i*10, stri);
outtextxy(ix,329,stri);
end;
ellipse (39, 320, 0, 90, 320, 320);
i := 0;
while i<10 do
begin
i := i+1;
fi := (i-1)*10*pi/180;
ix := round(scx*cos(fi)*30);
iy := 320-round(scy*sin(fi)*30);
jx := round((scx+0.4)*cos(fi)*30);
jy := 320-round((scy+0.3)*sin(fi)*30);
line(ix+40,iy,jx+40,jy);
str ((i-1)*10,stri);
outtextxy (40+jx+5,jy-5,stri);
end;
str(a[0,0]:10:3, stri);
outtextxy(350,35,' Radial ¯ '+stri);
str (a[30,0]:10:3, stri);
outtextxy (350, 55,' Tangential ¯ ' +stri);
str (a[30, 30]:10:3, stri);
outtextxy(350, 75,' Axial ¯ '+stri);
str(izotrop:10:3,stri);
outtextxy (350, 95,' Izotrop ¯ ' +stri);
outtextxy(0,340,'RADIAL');
outtextxy (0,10,'AXIAL' ) ;
outtextxy(330,340,'TANGENTIAL');
{ -----Plotarea liniilor de egal nivel---- }
settextstyle (0,0,2);
for m:=0 to m9 do
begin
nivel := minim+nivel0*m;
str (m:3,stri1);
stri1:='clasa '+stri1+' =';
str(nivel:10:6,stri);
outtextxy (420,120+20*m, stri1+stri);
end;
nivel:=minim;
for m:=1 to m9-1 do
begin
nivel:=minim+nivel0*m;
for i:=0 to 29 do
begin
for j:= 0 to 29 do
begin
n00:=a[i,j]-nivel;
nx1:=a[i+1,j]-nivel;
ny1:=a[i,j +1]-nivel;
n11:=a[i+1,j+1]-nivel;
if n00*nx1<0 then
{~intrare de sus~}
begin
transformare(i+n00/(n00-nx1),j,ix,iy);
moveto(ix,iy);
setcolor(1);
stanga;
jos;
dreapta;
setcolor(1);
end;
if n00*ny1<0 then
{~intrare din stanga~}
begin
transformare(i,j+n00/(n00-ny1),ix,iy);
moveto(ix,iy);
setcolor (1);
jos;
dreapta;
sus;
setcolor (1);
end;
if ny1*n11<0 then
{~intrare de jos}
begin
transformare(i+ny1/(ny1-n11),j+1,ix,iy);
moveto(ix,iy);
setcolor (1);
dreapta;
sus;
stanga;
end;
if nx1*n11<0 then
begin
transformare(i+1,j+nx1/(nx1-n11),ix,iy);
moveto(ix,iy);
setcolor(1);
sus;
jos;
stanga;
setcolor(1);
end;
end; { end j }
end;{ end i }
{readln;}{ ~sfarsitul plotarii unei clase~}
end;{ end m}
setcolor(1);
sound(1000) ;
delay(10);
nosound;
readln;
END.
|
Prof.Dr. D. Ciurchea
|
|
Linii de egal nivel-PASCAL
