6semforce (Вспомогательные материалы)

C Хайрулин И.Р. 06-421

C Содержит подпрограмму FORCE (упрощенная для дуги окружности)

C ПOДПPOГPAMMA PAЗHECEHИЯ HAГPУЗOK ПО УЗЛАМ

C BЫЗЫBAETCЯ ИЗ MAIN, вызываемых модулей нет.

C ======================================== ==========================

* IPR массив управления печатью

* CORD - одномерный массив глобальных координат узлов (первые два элемента этого массива

* являются координатами узла 1 по оси Х и Y соответственно,

* вторая пара элементов - узла 2 и т.д.).

* R - массив значений усилий в узлах (аналогичен по структуре массиву CORD,

* только вместо координат в нем указываются значения нагрузок по оси Х и Y в узлах);

* В примере:

* а) к дуге окружности радиусом 24 приложена равномерно распределенная нагрузка.

* б) к наклонной стороне с углом наклона 45 градусов приложена

* неравномерно распределенная нагрузка q(w) = q0*(sin(z*w))**3 ,

* имеющая положительное и отрицательное направление на разных участках,

*

C ======================================== ==========================

C wwwwwwwwwwwwwwwwww начало кода FORCE wwwwwwwwwwwwwwwwwwwwwwwwwwwww

C ======================================== ==========================

SUBROUTINE FORCE(NB,CORD,NRC,DB,IPR,NP,RSUM,R,NDF,D H,NR,

>PRM1,PRM2,PRM3,PRM4,PRM5,PRM6)

LOGICAL*1 IPR(50)

DIMENSION CORD(1),R(1),IFORCE1(100), IFORCE3(100),IFORCE4(100)

IF(IPR(22)) WRITE(6,22)PRM1,PRM4,PRM5,PRM2,PRM3

IF(IPR(23)) WRITE(6,20)

* I1,I3,I4 -кол-во точек к которым приложена расспределенная сила на участках

* NR - количество нагруженных узлов.

* NP - количество узлов

NR=0

I1=0

I3=0

I4=0

* обнуляем массивы

DO 123 I=1,100

IFORCE1(I)=0

IFORCE3(I)=0

123 IFORCE4(I)=0

*--------------------------------------- --------------------------

* !!Поиск узлов, к которым приложены силы!!

DO 100 I=1,NP

* Используем цикл по I до NP для обнуления массива сил

R(NDF*(I-1)+1)=0

R(NDF*(I-1)+2)=0

* Обрабатываемые точки

XT=CORD(NDF*(I-1)+1)

YT=CORD(NDF*(I-1)+2)

* Отлавливаем узлы на дуге радиусом 24 и записываем в массив IFORCE1

IF (ABS((XT-0)**2+(YT-60)**2-576).LT.думайт е) THEN

IF (ABS(PRM4-1).LT.0.01) THEN

I1=I1+1

* NR=NR+1

IFORCE1(I1)=I

ENDIF

ENDIF

* Отлавливаем узлы части наклонной, к которой приложена положительная нагрузка

* и записываем в массив IFORCE3

IF ((ABS(YT+XT-84).LT.0.1).AND.((XT-34.15). LT.0.1)) THEN

IF ((ABS(PRM6-2).LT.0.01).OR.(ABS(PRM6-1).L T.0.01)) THEN

I3=I3+1

* NR=NR+1

IFORCE3(I3)=I

ENDIF

ENDIF

* Отлавливаем узлы части наклонной, к которой приложена отрицательная нагрузка

* и записываем в массив IFORCE4

IF ((ABS(YT+XT-84).LT.0.1).AND.(((XT-34.15) .GT.0.1).OR.(ABS

> (XT-44.15).LT.0.01))) THEN

IF ((ABS(PRM6-3).LT.0.01).OR.(ABS(PRM6-1).L T.0.01)) THEN

I4=I4+1

* NR=NR+1

IFORCE4(I4)=I

ENDIF

ENDIF

100 CONTINUE

*--------------------------------------- ------------------------------

* !!Сортировки массивов с узлами!!

* Сортировка содержимого массива IFORCE1 по убыванию координаты Y

DO 201 I=I1,2,-1

DO 201 J=1,I-1

IF (CORD(NDF*(IFORCE1(J+1)-1)+2).GT.

>CORD(NDF*(IFORCE1(J)-1)+2)) THEN

ITMP=IFORCE1(J+1)

IFORCE1(J+1)=IFORCE1(J)

IFORCE1(J)=ITMP

201 ENDIF

* Сортировка содержимого массива IFORCE3 по убыванию координаты Y

DO 203 I=I3,2,-1

DO 203 J=1,I-1

IF (CORD(NDF*(IFORCE3(J+1)-1)+2).GT.

>CORD(NDF*(IFORCE3(J)-1)+2)) THEN

ITMP=IFORCE3(J+1)

IFORCE3(J+1)=IFORCE3(J)

IFORCE3(J)=ITMP

203 ENDIF

* Сортировка содержимого массива IFORCE4 по убыванию координаты Y

DO 204 I=I4,2,-1

DO 204 J=1,I-1

IF (CORD(NDF*(IFORCE4(J+1)-1)+2).GT.

>CORD(NDF*(IFORCE4(J)-1)+2)) THEN

ITMP=IFORCE4(J+1)

IFORCE4(J+1)=IFORCE4(J)

IFORCE4(J)=ITMP

204 ENDIF

*--------------------------------------- -----------------------------

* !!Задание сил на участках!!

* Сила приложенная к текущему узлу (вычисляется через площадь трапеции)

* S-размер участка поверхности между силами

* R1-распределенная сила приложенная к текущему узлу

* Задание сил в узлах по дуге окружности IFORCE1

DO 301 I=1,I1

J=IFORCE1(I)

* Координата текущей точки

X1=CORD(2*(IFORCE1(I)-1)+1)

Y1=CORD(2*(IFORCE1(I)-1)+2)

* Координата предыдущей точки

IF (I.NE.1) THEN

X3=CORD(2*(IFORCE1(I-1)-1)+1)

Y3=CORD(2*(IFORCE1(I-1)-1)+2)

ELSE

X3=X1

Y3=Y1

ENDIF

* Координата следующей точки

IF (I.NE.I1) THEN

X2=CORD(2*(IFORCE1(I+1)-1)+1)

Y2=CORD(2*(IFORCE1(I+1)-1)+2)

ELSE

X2=X1

Y2=Y1

ENDIF

* Высчитываем S - длину основания, на которую будет действовать нагрузка,

* приходящаяся на узел.

* В формуле трапеций S - высота трапеции,

* для окружностей - длина дуги (аппроксимированной.Равна длине хорды для упрощения)

* подпрограмма- функция OSNOV подсчитывает основание как сумму половин хорд

* (расстояний) от данного узла до соседних узлов.

S=OSNOV(X1,Y1,X2,Y2,X3,Y3)

* Рассчитываем R1 значение интеграла (силы, приложенной к узлу).

* PRM1- значение равномерно распределенной нагрузки по дуге окружности.

*

R1=S*PRM1

* Заносим в массив R получившиеся нагрузки разложенные по X и Y

R(NDF*(IFORCE1(I)-1)+1)=R(NDF*(IFORCE1(I )-1)+1)-R1*(X1-0)/

>sqrt(ABS((X1-0)**2+(Y1-60)**2))

R(NDF*(IFORCE1(I)-1)+2)=R(NDF*(IFORCE1(I )-1)+2)-R1*(Y1-60)/

>sqrt(ABS((X1-0)**2+(Y1-60)**2))

* Выводим нагруженные узлы, их координаты и подсчитываем кол-во

IF(IPR(23)) THEN

IF((ABS(PRM6-1).GT.0.01).AND.(ABS(PRM6-2 ).GT.0.01)) THEN

WRITE(6,21) J,CORD(NDF*(J-1)+1),

> CORD(NDF*(J-1)+2),R((J-1)*NDF+1),R((J-1) *NDF+2)

NR=NR+1

ELSE IF ((ABS(CORD(NDF*(J-1)+1)-24).GT.0.1).AND.

> (ABS(CORD(NDF*(J-1)+2)-60).GT.0.1)) THEN

WRITE(6,21) J,CORD(NDF*(J-1)+1),

> CORD(NDF*(J-1)+2),R((J-1)*NDF+1),R((J-1) *NDF+2)

NR=NR+1

ENDIF

ENDIF

301 continue

* Задание сил по узлам от неравномерно распределенной нагрузки

* на участок, где нагрузка приложена в положительном направлении

* IFORCE3 часть наклонной, к которой приложена положительная нагрузка

DO 303 I=1,I3