12PosobKR2015 (1013897), страница 11
Текст из файла (страница 11)
R(NDF*(IFORCE2(I)-1)+2)=R(NDF*(IFORCE2(I)-1)+2)
>+R1*SQRT(1/(1+(3.)**2))
ELSE
! для переноса нагрузки находим коэф(нужно только для нижней точки)
R2=(R1-R1*SQRT(1-1/(1+(3.)**2)))
! перенос нагрузки:
R(NDF*(IFORCE2(I)-1)+1)=R(NDF*(IFORCE2(I)-1)+1)
>+R1*SQRT(1-1/(1+(3.)**2)) ! коэф. 3 взят из тангенса(коэф. - k) уравнения прямой
R(NDF*(J_PRED)+2)=R(NDF*(J_PRED)+2)+R2*SQRT(1/(1+(3.)**2))
R(NDF*(IFORCE2(I)-1)+1)=R(NDF*(IFORCE2(I)-1)+1)
>+R2*SQRT(1/(1+(3.)**2))
ENDIF
! Вывод текстовой информации
IF(IPR(23)) WRITE(6,21) J,CORD(NDF*(J-1)+1),
> CORD(NDF*(J-1)+2),R((J-1)*NDF+1),R((J-1)*NDF+2)
400 CONTINUE
IF(IPR(23)) PRINT *,"Число нагруженных узлов NR=",NR
21 FORMAT (I6,1X,2F10.3,4X,2F9.2)
122 RETURN
IF (NR.NE.0) GOTO 122
WRITE(6,6) (I,(R((I-1)*NDF+J),J=1,NDF),I=1,NP)
6 FORMAT (' Не найдены точки приложения нагрузок'/8(I4,2G10.2))
37 FORMAT (F10.3)
20 FORMAT (' Распределение сил/ ' Узел X', ' Y * RX * * RY *')
22 FORMAT (' Суммарная сила '/' * RX * * RY *'/ (3X,2F11.3))
STOP
END
*Подпрограмма вычисления интеграла методом трапеций для прямой Y=180-3X
FUNCTION RPLOSH(XP,XT,XN,PAR) !входные значения(координаты X) - предыдущая(XP), текущая(XT), следующая(XN)
RPLOSH=0
F_COS = SQRT(1-1/(1+(3.)**2))
XFN = ( ((60.-XT)+(60.-XP))/2.)*3./F_COS
XFK = ( ((60.-XT)+(60.-XN))/2.)*3./F_COS
SH=(XFK-XFN)/10.
XT=XFN+SH
DO 25 IS=1,10,1
RPLOSH=RPLOSH+SH*( (((PAR-(XT-SH))**1.5)*COS(0.05*(XT-SH)))+
>(((PAR-XT)**1.5)*COS(0.05*XT)))/2.
XT=XT+SH
25 CONTINUE
RETURN
END
* ==========================
* конец кода FORCE
*===========================
Пояснение к работе RPLOSH
На вход RPLOSH (X3,X1,X2, PRM3) подаются:
X1 - глобальная координата узла по оси х, к которому приложена сила, подсчитываемая в RPLOSH через интеграл от функции нагрузки q(Пси) = PAR –((Пси)**1.5)*Cos(0.05*Пси);
Х3 - глобальная координата Х предыдущего узла;
X2 глобальная координата Х следующего узла.
XFN и XFK в RPLOSH – это точки начала и конца интегрирования, записанные через локальную координату Пси (отсчитываемую от начала прямой стороны - её нижней точки), по которой действует функция нагрузки q(Пси). Эти XFN и XFK находятся на середине расстояния между локальными координатами текущего и предыдущего узла и середине расстояния между локальными координатами текущего и следующего узла.
Так как координата Пси=0 находится в нижней точке наклонной стороны (см. рис. 43), имеющей глобальные координаты (60:0), то XFN и XFK подсчитываются с учетом значения 60.
В данной подпрограмме отрезок [XFN ;XFK] делится на 10 частей, после чего вычисляется значение локальной координаты XT, которую не надо путать с поступившем в RPLOSH значении глобальной координаты XT=Х1. Просто разработчик, использовал уже не нужный ему идентификатор XT второй раз, хотя мог бы поступить иначе.
Результат работы подпрограммы Force в подсистеме графического отображения результатов расчета будет представлен в следующем виде (рис. 46)
Обратите внимание, что в самом нижнем правом узле сила не перпендикулярна стороне, к которой она приложена. Это связано с тем, что в случае правильного приложения силы она бы имела составляющую вдоль оси Y. Но рассматриваемый узел закреплен от перемещения вдоль оси Y. МКЭ не допускает приложения силы к закрепленному узлу в направлении оси, по которой запрещено перемещение (программа просто не будет правильно считать). Поэтому здесь приходится прибегать к искусственному приему, исключая вертикальную составляющую действующей на этот узел силы. Ещё более часто всю силу переносят на соседний незакрепленный узел. При увеличении числа конечных элементов образующаяся ошибка уменьшается.
В текстовом файле результатов информация о нагруженных узлах будет выведена в следующем виде:
*****Start FORCE (реализация внешних воздействий)
2 .000 60.000 .00 11327.08
7 3.545 60.000 .00 20525.93
12 7.090 60.000 .00 17910.79
17 10.635 60.000 .00 15562.96
21 14.180 60.000 .00 12719.04
22 17.295 60.000 .00 10432.22
23 20.410 60.000 .00 9201.02
24 23.525 60.000 .00 8151.18
25 26.640 60.000 .00 7531.77
26 29.980 60.000 .00 7027.51
27 33.320 60.000 .00 6468.61
38 36.660 60.000 .00 6133.28
137 60.000 .000 9883.97 .00
136 58.606 4.190 19797.02 6765.23
134 57.030 8.920 19717.79 6572.60
131 55.271 14.190 17378.81 5792.94
126 53.330 20.000 11260.22 3753.41
119 51.694 24.906 5048.45 1682.82
110 50.060 29.810 -1.39 -.46
100 48.429 34.711 -4724.57 -1574.86
89 46.800 39.610 -9068.57 -3022.86
77 45.100 44.711 -12763.16 -4254.39
68 43.400 49.810 -15161.07 -5053.69
59 41.700 54.906 -16313.91 -5437.97
49 40.000 60.000 -8198.11 277.98
Число нагруженных узлов NR= 25
*****Finish FORCE
2.7. Отладка расчетного блока и получение результатов
После составления подпрограмм Bound и Force, занесения информации о свойствах КЭ и вспомогательных значений в окно основных параметров (рис.45) проект готов к компиляции. Если компиляция прошла успешно, то можно произвести расчет напряжений.
Результаты расчета можно просмотреть как в численном, так и в графическом виде. Для просмотра численных результатов расчета нажмите пиктограмму «Текстовый результат» или найдите такую альтернативу в выпадающем меню «Результат». Для просмотра результатов в графическом виде нажмите кнопку «Графическое отображение результатов расчета» или тоже воспользуйтесь выпадающим меню.
Графическое изображение результатов расчета может быть цветным (рис. 46) или чернобелым (рис. 47). В подсистеме графического вывода имеются настройки, управляющие цветом.
Проверьте правильность задания исходных данных с помощью анализа численных результатов расчета и графического изображения:
-
рисунок пластины в графическом виде должен напоминать исходную пластину;
-
сверьте координаты закрепленных узлов, а также тип закрепления;
-
проверьте по координатам, чтобы нагрузка была приложена к нужным узлам, а в сумме нагрузка соответствовала суммарной.
Если вы не увидели картинки напряжений, то, очевидно, вы сделали алгоритмическую ошибку в модифицированных подпрограммах или исходных данных. Для определения этой ошибки просмотрите окно с численными результатами расчета. Если информации в нем недостаточно, то просмотрите любым способом полный файл результатов расчета, имя которого вы задавали как имя файла результатов расчета. Он хранится в каталоге SIGMA\DEBUG. Определив эти ошибки и исправив их, откомпилируйте проект и попробуйте снова рассчитать напряжения.
Попытайтесь найти конечный элемент содержащий какую-либо определенную точку, а также напряжение в этом элементе. Для этого в численных результатах найдите номер конечного элемента, содержащего выбранную точку. Сделать это будет легче, если вы будете ориентироваться на расположенный рядом узел зоны. Затем по номеру элемента найдите в таблице значение его напряжения.
Чтобы ваши упражнения несли прикладной характер, найдите элемент с наибольшим напряжением и, варьируя либо толщину пластины, либо нагрузку, добейтесь разрушения пластины – когда напряжение в элементе превысит предельное напряжение для материала, из которого выполнена пластина.
Список литературы
-
А.Джордж, Дж. Лю., Численное решение больших разреженных систем уравнений. - М.: МИР, 1984.
-
Jim Ruppert: A Delaunay Refinement Algorithm for Quality 2-Dimentional Mesh Generation, NASA Ames Research Center, Submission to Journal of Algorithms, 1994.
-
Л. Сегерлинд., Применение метода конечных элементов. - М.: МИР, 1979.