4. Zastosowanie MES do modelowania odkształcenia sprężysto-plastycznego

1. Podstawy teoretyczne

Rys. Idea metody sprężystych rozwiązani

(4.1)

2. Dyskretyzacja ośrodka

Rys.2. Przykład siatki MES

Rys. 3. Rozpatrywany przykład siatki MES

Rys. 4. Numeracja węzłów i elementów

Funkcjonał po dyskretyzacji

3. Struktury danych

ne=6

nh=7

Ux(1)=0; Uy(1)=Uup;

Ux(2)=0; Uy(2)=Uup;

Ux(3)=0; Uy(3)=0;

Ux(4)=0; Uy(4)=Udown;

Ux(5)=0; Uy(5)=Udown;

Ux(6)=0; Uy(6)=0;

Ux(7)=0; Uy(7)=0 ;

Status1 = 1

Status2 = 1

Status3 = 0

Status4 = 1

Status5 = 1

Status6 = 0

Status7 = 0

X1=-1 Y1=1

X2=1 Y2=1

X3=1 Y3=0

X4=1 Y4=-1

X5=-1 Y5=-1

X6=-1 Y6=0

X7=0 Y7=0

nop(1,1)=1; nop(2,1)=7; nop(3,1)=2;

nop(1,2)=7; nop(2,2)=3; nop(3,2)=2;

nop(1,3)=3; nop(2,3)=7; nop(3,3)=4;

nop(1,4)=7; nop(2,4)=5; nop(3,4)=4;

nop(1,5)=7; nop(2,5)=6; nop(3,5)=5;

nop(1,6)=6; nop(2,6)=7; nop(3,6)=1.

4. Algorytm obliczenia funkcjonału

5. Program

'***************************************************

' Program do modelowania spęczania

' za pomocą MES

' Milenin A.

' milenin@agh.edu.pl

'***************************************************

' Constants

Const g_nh = 500, g_ne = 1000

Type node

x As Double

y As Double

Ux As Double

Uy As Double

Ex As Double

Ey As Double

Exy As Double

Sx As Double

Sy As Double

Sxy As Double

Ei As Double

Si As Double

status As Integer

LOK(0 To 8) As Integer

End Type

Type element

nop(1 To 3) As Integer

End Type

Type dataFEM

nh As Integer

ne As Integer

dVel As Double

Gs As Double

k As Double

Uup As Double

Udown As Double

H0 As Double

B0 As Double

nhB As Integer

nhH As Integer

ND(1 To g_nh) As node

EL(1 To g_ne) As element

End Type

Dim gr As dataFEM

'***********************************************

' Program glowny

'***********************************************

Sub mainFEM()

Call DataIni

' Call DataGen

Call DataGenAuto

Call GetLOK

Call VisualData(0)

Call solve

Call VisualData(1)

End Sub

'*************************************************

Sub solve()

Dim ih As Integer

Dim iter As Integer

Dim F0 As Double

Dim Fplus As Double

Dim dF As Double

For iter = 1 To 800

For ih = 1 To gr.nh

If (gr.ND(ih).status = 0) Then

F0 = GetFunNode(ih)

gr.ND(ih).Ux = gr.ND(ih).Ux + gr.dVel

Fplus = GetFunNode(ih)

If (Fplus > F0) Then

gr.ND(ih).Ux = gr.ND(ih).Ux - 2 * gr.dVel

End If

gr.ND(ih).Uy = gr.ND(ih).Uy + gr.dVel

Fplus = GetFunNode(ih)

If (Fplus > F0) Then

gr.ND(ih).Uy = gr.ND(ih).Uy - 2 * gr.dVel

End If

If (Abs(F0) < 0.000000001) Then F0 = 0.000000001

dF = Abs((Fplus - F0) / F0)

End If

Next ih

Next iter

For ih = 1 To gr.nh

gr.ND(ih).x = gr.ND(ih).x + gr.ND(ih).Ux

gr.ND(ih).y = gr.ND(ih).y + gr.ND(ih).Uy

Next ih

End Sub

Function GetFunNode(ih As Integer)

Dim i As Integer

Dim ie As Integer

Dim Je As Double

GetFunNode = 0

For i = 1 To gr.ND(ih).LOK(0)

ie = gr.ND(ih).LOK(i)

Je = clcFunkconal(ie)

GetFunNode = GetFunNode + Je

Next i

End Function

Sub GetLOK()

Dim ih As Integer

Dim ie As Integer

Dim j As Integer

For ih = 1 To gr.nh

gr.ND(ih).LOK(0) = 0

For ie = 1 To gr.ne

For j = 1 To 3

id = gr.EL(ie).nop(j)

If (id = ih) Then

gr.ND(ih).LOK(0) = gr.ND(ih).LOK(0) + 1

gr.ND(ih).LOK(gr.ND(ih).LOK(0)) = ie

End If

Next j

Next ie

Next ih

End Sub

Function clcFunkconal(ie As Integer)

Dim x(1 To 3) As Double

Dim y(1 To 3) As Double

Dim Ux(1 To 3) As Double

Dim Uy(1 To 3) As Double

Dim Ex As Double

Dim Ey As Double

Dim Exy As Double

Dim Ei As Double

Dim E0 As Double

Dim Si As Double

Dim E As Double

Dim k As Double

Dim Sx As Double

Dim Sy As Double

Dim Sxy As Double

Dim i(1 To 3) As Integer

Dim bi As Double

Dim bj As Double

Dim bk As Double

Dim ci As Double

Dim cj As Double

Dim ck As Double

Dim Ae As Double

Dim a As Double

Dim b As Double

Dim Je As Double

x(1) = gr.ND(gr.EL(ie).nop(1)).x

x(2) = gr.ND(gr.EL(ie).nop(2)).x

x(3) = gr.ND(gr.EL(ie).nop(3)).x

y(1) = gr.ND(gr.EL(ie).nop(1)).y

y(2) = gr.ND(gr.EL(ie).nop(2)).y

y(3) = gr.ND(gr.EL(ie).nop(3)).y

Ux(1) = gr.ND(gr.EL(ie).nop(1)).Ux

Ux(2) = gr.ND(gr.EL(ie).nop(2)).Ux

Ux(3) = gr.ND(gr.EL(ie).nop(3)).Ux

Uy(1) = gr.ND(gr.EL(ie).nop(1)).Uy

Uy(2) = gr.ND(gr.EL(ie).nop(2)).Uy

Uy(3) = gr.ND(gr.EL(ie).nop(3)).Uy

bi = y(2) - y(3)

ci = x(3) - x(2)

bj = y(3) - y(1)

cj = x(1) - x(3)

bk = y(1) - y(2)

ck = x(2) - x(1)

Ae = x(2) * y(3) + x(1) * y(2) + y(1) * x(3) - y(1) * x(2) - y(2) * x(3) - x(1) * y(3)

Ex = (bi * Ux(1) + bj * Ux(2) + bk * Ux(3)) / (2 * Ae)

Ey = (ci * Uy(1) + cj * Uy(2) + ck * Uy(3)) / (2 * Ae)

Exy = (bi * Uy(1) + bj * Uy(2) + bk * Uy(3) + ci * Ux(1) + cj * Ux(2) + ck * Ux(3)) / (4 * Ae)

Ei = (((2) ^ 0.5) / 3) * ((Ex - Ey) ^ 2 + Ex ^ 2 + Ey ^ 2 + 6 * Exy ^ 2)

E0 = (Ex + Ey) / 3

' Sx

' Sy

' Sz

a = 2000

b = 0.18

Si = a * Ei ^ b ' krzywa umocnienia

E = 1 ' debug

k = 1000 ' debug

' Je = a * Ei ^ (b + 1) / (b + 1) * Ae + k * E0 ^ 2 * Ae

Je = E * Ei ^ 2 * Ae + k * E0 ^ 2 * Ae

clcFunkconal = Je

End Function

Sub DataIni()

Dim i As Integer

For i = 1 To g_nh

gr.ND(i).x = 0

gr.ND(i).y = 0

gr.ND(i).x = 0

gr.ND(i).y = 0

gr.ND(i).Sy = 0

gr.ND(i).Sxy = 0

gr.ND(i).Sx = 0

gr.ND(i).Si = 0

gr.ND(i).Ei = 0

gr.ND(i).Ex = 0

gr.ND(i).Exy = 0

gr.ND(i).Ey = 0

gr.ND(i).status = 0

Next i

End Sub

Sub DataGenAuto()

Dim x As Double

Dim y As Double

Dim dx As Double

Dim dy As Double

Dim i As Integer

Dim j As Integer

Dim id As Integer

Dim n1 As Integer

Dim n2 As Integer

Dim n3 As Integer

Dim n4 As Integer

' Dane podstawowe

gr.dVel = 0.001

gr.Uup = -1

gr.Udown = 1

gr.B0 = 10 'mm

gr.H0 = 10 'mm

gr.nhB = 6

gr.nhH = 6

gr.nh = gr.nhB * gr.nhH

gr.ne = (gr.nhB - 1) * (gr.nhH - 1) * 2

dx = gr.B0 / (gr.nhB - 1)

dy = gr.H0 / (gr.nhH - 1)

y = 0

id = 1

For j = 1 To gr.nhH

x = 0

For i = 1 To gr.nhB

gr.ND(id).x = x

gr.ND(id).y = y

gr.ND(id).status = 0

If (j = 1) Then

gr.ND(id).status = 1

gr.ND(id).Ux = 0

gr.ND(id).Uy = gr.Udown

End If

If (j = gr.nhH) Then

gr.ND(id).status = 1

gr.ND(id).Ux = 0

gr.ND(id).Uy = gr.Uup

End If

x = x + dx

id = id + 1

Next i

y = y + dy

Next j

id = 1

For j = 1 To gr.nhH - 1

For i = 1 To gr.nhB - 1

ih = j + (i - 1) * (gr.nhH - 1)

'

' n4 *---------* n1

' I !

' I !

' I !

' I !

' n3 *---------* n2

'

n1 = j + (i - 1) * gr.nhH + gr.nhH + 1

n2 = j + (i - 1) * gr.nhH + 1

n3 = j + (i - 1) * gr.nhH

n4 = j + (i - 1) * gr.nhH + gr.nhH

gr.EL(id).nop(1) = n1

gr.EL(id).nop(2) = n3

gr.EL(id).nop(3) = n2

gr.EL(id + 1).nop(1) = n1

gr.EL(id + 1).nop(2) = n4

gr.EL(id + 1).nop(3) = n3

id = id + 2

Next i

Next j

End Sub

Sub DataGen()

gr.dVel = 0.001

gr.Gs = 1

gr.k = 100

gr.Uup = -0.1

gr.Udown = 0.1

gr.nh = 7

gr.ne = 6

gr.ND(1).x = -10

gr.ND(1).y = 10

gr.ND(1).status = 1

gr.ND(1).Ux = 0

gr.ND(1).Uy = gr.Uup

gr.ND(2).x = 10

gr.ND(2).y = 10

gr.ND(2).status = 1

gr.ND(2).Ux = 0

gr.ND(2).Uy = gr.Uup

gr.ND(3).x = 10

gr.ND(3).y = 0

gr.ND(3).status = 0

gr.ND(3).Ux = 0

gr.ND(3).Uy = 0

gr.ND(4).x = 10

gr.ND(4).y = -10

gr.ND(4).status = 1

gr.ND(4).Ux = 0

gr.ND(4).Uy = gr.Udown

gr.ND(5).x = -10

gr.ND(5).y = -10

gr.ND(5).status = 1

gr.ND(5).Ux = 0

gr.ND(5).Uy = gr.Udown

gr.ND(6).x = -10

gr.ND(6).y = 0

gr.ND(6).status = 0

gr.ND(6).Ux = 0

gr.ND(6).Uy = 0

gr.ND(7).x = 0

gr.ND(7).y = 0

gr.ND(7).status = 0

gr.ND(7).Ux = 0

gr.ND(7).Uy = 0

gr.EL(1).nop(1) = 1

gr.EL(1).nop(2) = 7

gr.EL(1).nop(3) = 2

gr.EL(2).nop(1) = 7

gr.EL(2).nop(2) = 3

gr.EL(2).nop(3) = 2

gr.EL(3).nop(1) = 3

gr.EL(3).nop(2) = 7

gr.EL(3).nop(3) = 4

gr.EL(4).nop(1) = 7

gr.EL(4).nop(2) = 5

gr.EL(4).nop(3) = 4

gr.EL(5).nop(1) = 7

gr.EL(5).nop(2) = 6

gr.EL(5).nop(3) = 5

gr.EL(6).nop(1) = 6

gr.EL(6).nop(2) = 7

gr.EL(6).nop(3) = 1

End Sub

Sub VisualData(kod As Integer)

Dim i As Integer

Dim pnt1(0 To 2) As Double

Dim pnt2(0 To 2) As Double

Dim pnt3(0 To 2) As Double

Dim color As AcadAcCmColor

Dim myGrid As Acad3DFace

For i = 1 To gr.ne Step 1

pnt1(0) = gr.ND(gr.EL(i).nop(1)).x

pnt1(1) = gr.ND(gr.EL(i).nop(1)).y

pnt1(2) = 0

pnt2(0) = gr.ND(gr.EL(i).nop(2)).x

pnt2(1) = gr.ND(gr.EL(i).nop(2)).y

pnt2(2) = 0

pnt3(0) = gr.ND(gr.EL(i).nop(3)).x

pnt3(1) = gr.ND(gr.EL(i).nop(3)).y

pnt3(2) = 0

Set myGrid = ThisDrawing.ModelSpace.Add3DFace(pnt1, pnt2, pnt3, pnt3)

Set color = AcadApplication.GetInterfaceObject("AutoCAD.AcCmColor.16")

If (kod = 0) Then

Call color.SetRGB(80, 100, 244)

Else

Call color.SetRGB(128, 100, 128)

End If

myGrid.TrueColor = color

Next i

End Sub

6. Wyniki

6. Problemy

7. Sprawozdanie

Tekst programu w dowolnym języku oprogramowania z obliczonym przykładem testowym

Opis algorytmów programu w tekście (komentarzy)