| Правила | Регистрация | Пользователи | Поиск | Сообщения за день | Все разделы прочитаны |  Справка по форуму | Файлообменник |

Вернуться   Форум DWG.RU > Программное обеспечение > Программирование > Как создать расчетное программное обеспечение с открытым исходным кодом (конструктивные решения)

Как создать расчетное программное обеспечение с открытым исходным кодом (конструктивные решения)

Ответ
Поиск в этой теме
Непрочитано 24.09.2021, 14:52
Как создать расчетное программное обеспечение с открытым исходным кодом (конструктивные решения)
nickname2019
 
Регистрация: 18.11.2019
Сообщений: 1,039

На мой взгляд, основными проблемами российского рынка расчетного программного обеспечения являются:
- отсутствие нормальной возможности программной автоматизации по решению расчетных задач
(в расчетных программах отсутствует возможность для нормального программирования, т.е. невозможно написать программу для полностью автоматического создания расчетной схемы (нескольких расчетных схем), автоматического выполнения расчета, автоматического получения результатов и их автоматического анализа);
- закрытый исходный код по подбору расчетных параметров несущих элементов
(различные программы дают различные результаты при решении одинаковых задач, сравнение алгоритмов подбора различных между собой невозможно, так как код закрыт, общепринятых и одобренных алгоритмов нет, каждый пользуется своим "черным ящиком", который иногда может выдать ошибочное решение);
- для людей, которые занимаются автоматизацией на Лисп, C# и т.д. отсутствуют инструменты, которые позволяли бы программно "по-простому" вызвать готовую библиотечную функцию (например, по подбору сечения какой-то простой балки непосредственно из графического редактора), что вызывает необходимость вызова отдельной расчетной программы, что серьезно тормозит работу;
- "корявый" интерфейс, ужасно неудобная и медленная работа в существующих российских (и украинских) расчетных программах;
(фактически при наличии нормального графического редактора (autocad, nanocad и т.д.) приходится экспортировать данные в в "корявый" редактор расчетной программы и длительное в нем работать (задавать нагрузки, связи и т.д.), а встроить расчетную программу в нормальный графический редактор через автоматизацию невозможно).

В связи с вышеизложенным, назрел вопрос:
Как технологически наиболее правильно можно организовать разработку расчетного программного обеспечения с открытым исходным кодом?

Для совместной разработки кода создано общее хранилище на GitHub, используя которое каждый может поучаствовать в разработке :
https://github.com/chaosEagleOwl/source

На данным момент работа находится в стадии тестирования возможности совместной разработки.
Требования к программному обеспечению изложены в файле (ссылка README.md на GitHub): https://github.com/chaosEagleOwl/source/README.md

ТЗ на модуль формирования КЭ-сеток сформировано и помещено на GitHub.
На весь комплекс ТЗ формировать долго, видимо, будет чуть позже.

Сформирована доска для управления проектом, туда добавлены наиболее актуальные задачи.
Задачи проекта.

Последний раз редактировалось nickname2019, 06.10.2021 в 09:07.
Просмотров: 42962
 
Непрочитано 24.10.2021, 13:25
#321
bigden


 
Регистрация: 05.08.2008
Сообщений: 801


я конечно не спец в линейной алгебре, но

Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Матрицы держит в памяти (вроде даже без учета симметрии - хранит и верх, и низ). Значит, большие схемы на слабом железе считать не выйдет?
держать полностью матрицу по времени должно быстрее получаться, чем прыгать по треугольнику и на каждой операции выгружать и загружать в кеш. в каждой древней книге написано, что нужно хранить только половину. я подозреваю, что это пережиток прошлого, когда компьютеры были большие, а память маленькая.

Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Не понимает, что ноль на диагонали - это автоматическая связь, а не ошибка расчета. Можно ему объяснить, что он неправ?
что это за фишка такая - расчет с нулем на диагонали? можно где-то посмотреть?
bigden вне форума  
 
Непрочитано 24.10.2021, 13:39
#322
zvezdochiot

маркшейдер
 
Регистрация: 25.09.2021
Москва
Сообщений: 37


Цитата:
Сообщение от bigden Посмотреть сообщение
я подозреваю, что это пережиток прошлого, когда компьютеры были большие, а память маленькая.
Да не только по алгебре неспец. Вообще в памяти не хранились, а каждый раз считывались с носителя. Так что мимо!
__________________
Keep it simple, stupid.
zvezdochiot вне форума  
 
Непрочитано 25.10.2021, 04:42
1 | 1 #323
Нубий-IV

Инженер-философ
 
Регистрация: 24.04.2019
Хабаровск
Сообщений: 1,041


Цитата:
Сообщение от bigden Посмотреть сообщение
держать полностью матрицу по времени должно быстрее получаться
Так ясен пень, что быстрее, если памяти хватает.

Простой тест: типовая свечка на двадцать этажей, схема в 120000 узлов - это совсем скромно; старк считает ее две минуты, из них решение системы - около 10сек; остальное время занимают чтение данных и запись результатов.
В протоколе указано число ненулевых элементов матрицы - 0.238 млрд. Если считать, что матрица хранится в виде значений double (8 байт) с двумя индексами int (2x4 байта), это требует расхода памяти 0.238*(8+2*4) = 3.8Gb. Примерно столько и лежит у старка в рабочем каталоге; пиковое потребление памяти в процессе счета у него меньше, т.е. матрицу он выгружает на диск.

Тот же тест в Eigen.
Проверил: при решении через SimplicialLDLT не нужно записывать нижний треугольник, т.е. избыточное хранение есть только в плотных матрицах - у нас такими будут разве что матрицы суперэлементов, если мы их оставим с типом MatrixXd, как сейчас сделано в Types.h.
В верхнем треугольнике матрицы размером N, с шириной ленты W, и оставшейся нулевой шириной Z=N-W хранится элементов. Значит, ширина ленты равна
Код:
[Выделить все]
 
#include <iostream>
#include <ctime>

#include <stdio.h>
#include <tchar.h>
#include <math.h>

#include <Eigen/Dense>
#include <Eigen/Sparse>

int _tmain(int argc, _TCHAR* argv[])
{
	int MatrixSize = 716229;
	int BandWidth = 333;

	std::cout << "Stiffness matrix calulation...\n";
	unsigned start_time = std::clock();

	Eigen::SparseMatrix<double> K(MatrixSize, MatrixSize);
	std::vector<Eigen::Triplet<double>> T;
	for(int i = 0; i < MatrixSize; i++) {
		int bandMinIndex = i;
		int bandMaxIndex = std::min(i + BandWidth, MatrixSize);
		for(int j = bandMinIndex; j < bandMaxIndex; j++) {
			double distance = i-j;
			double max_distance = BandWidth;
			double relativeDistance = distance / max_distance;
			T.push_back(Eigen::Triplet<double>(i, j, 1.0 - std::pow(relativeDistance, 2)));
		}
	}
	K.setFromTriplets(T.begin(), T.end());

	unsigned end_time = std::clock();
	std::cout << (end_time - start_time) / 1000.0 << " seconds\n";
	std::cout << "...done.\n\n";
	//K.coeffRef(0,0) = 0.0;

	std::cout << "Forces vector calulation...\n";
	start_time = std::clock();
	Eigen::VectorXd F(MatrixSize);

	for(int i = 0; i < MatrixSize; i++)
		F[i] = (double)i / MatrixSize;

	end_time = std::clock();
	std::cout << (end_time - start_time) / 1000.0 << " seconds\n";
	std::cout << "...done.\n\n";

	std::cout << "Equation system solving...\n";
	start_time = std::clock();

	Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>> solver(K);
	Eigen::VectorXd x = solver.solve(F);

	end_time = std::clock();
	std::cout << (end_time - start_time) / 1000.0 << " seconds\n";
	std::cout << "...done.\n";

	//std::cout << "Solver info: " << solver.info() << std::endl;
	//std::cout << "\nStiffness:\n" << K << std::endl;
	//std::cout << "\nForces:\n" << F << std::endl;
	//std::cout << "\nDisplacements:\n" << x << std::endl;

	return 0;
}
Если собрать 32-битный проект, результат запуска будет таким:
Старенькие автокады x32 отпадают. На x64 результат получше:
Потрачено 7.4G памяти; у меня на рабочем компе память 16G, ее хватает. Забавно, что сборка матрицы шла в 10 раз дольше, чем решение системы; загрузка процессора 13% - работает 1 ядро из 8; и это без вычисления матриц жесткости элементов, просто запись чисел. То ли сборку надо распараллеливать, то ли операция добавления в матрицу тормозит; надо тестировать дополнительно.

Чтобы увеличить расход памяти, удвою ленту: W = 2*333 = 666. Результат:
После 80% загрузки памяти начинается использование подкачки. Память, похоже, запрашивается блоками. Каждый раз при запросе (зубцы на графике) на десяток-другой секунд начинаются тормоза типа "мышь повесилась". Однако, расчет выполняется; сборка матрицы втрое дольше, решение - вдвое; это вполне приличный результат.
Цитата:
Сообщение от bigden Посмотреть сообщение
я конечно не спец в линейной алгебре
Это не страшно. Через десять лет, появится наконец первая прилично работающая версия, все участники проекта станут спецами: в процессе мы просто обязаны собрать все грабли, предусмотренные этой наукой за прошедший век .
Миниатюры
Нажмите на изображение для увеличения
Название: 01.png
Просмотров: 286
Размер:	90.8 Кб
ID:	242040  Нажмите на изображение для увеличения
Название: 02.png
Просмотров: 294
Размер:	50.6 Кб
ID:	242041  Нажмите на изображение для увеличения
Название: 03.png
Просмотров: 295
Размер:	79.6 Кб
ID:	242042  Нажмите на изображение для увеличения
Название: 04.png
Просмотров: 291
Размер:	79.4 Кб
ID:	242043  
Нубий-IV вне форума  
 
Автор темы   Непрочитано 25.10.2021, 08:25
#324
nickname2019


 
Регистрация: 18.11.2019
Сообщений: 1,039


Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Не понимает, что ноль на диагонали - это автоматическая связь, а не ошибка расчета. Можно ему объяснить, что он неправ?
При формировании глобальной матрицы не следует нумеровать степени свободы, на которые наложены связи. Т.е. в глобальную матрицу жесткости эти степени свободы не должны попадать.
Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Приостановка/возобновление расчета невозможны?
Нужно в цикле проверять нажатие Esc. В автокаде есть примерно такой код.
Код:
[Выделить все]
 
acutPrintf(ACRX_T("\nФормирование матрицы...\r"));
 
acedSetStatusBarProgressMeter(ACRX_T("Формирвоание матрицы"), 0, n);

for (i=0;i<n;i++){

  if ( acedUsrBrk() ){
	// Perform cleanup on cancel
	bContinue = Adesk::kFalse;
	acutPrintf(ACRX_T("\nCancelled !"));
	break;
  }
  acedSetStatusBarProgressMeterPos(i);
//Делаем что-то
//....

}
acedRestoreStatusBar();

В посте #320 матрицу жесткости какого КЭ вычисляли? Этот код можно забрать в проект?


P.S. Судя по затратам времени на вычисления, делать класс - прокладку между eigen и основной программой - не самая удачная идея. Прокладка, видимо, будет тормозить расчет, так как количество вызовов процедур сразу в два раза вырастет (выделение/освобождение стека и т.д.).

При вычислении локальных матриц жесткостей, видимо, лучше использовать полностью-заполненные матрицы, при этом выделять память под локальную матрицу один раз, а заполнять ее столько раз, сколько есть КЭ данного типа в схеме.

Последний раз редактировалось nickname2019, 25.10.2021 в 08:55.
nickname2019 вне форума  
 
Непрочитано 25.10.2021, 08:58
#325
румата


 
Регистрация: 06.04.2015
Сообщений: 1,515


Цитата:
Сообщение от nickname2019 Посмотреть сообщение
При вычислении локальных матриц жесткостей, видимо, лучше использовать полностью-заполненные матрицы, при этом выделять память под локальную матрицу один раз, а заполнять ее столько раз, сколько есть КЭ данного типа в схеме.
Плотные матрицы, в принципе, нужно использовать там где это возможно. Скорость работы с ними значительно выше, чем с разреженными. Разреженные матрицы есть смысл использовать тогда когда экономия памяти становится значительно важней скорости работы, т.е. тогда, когда размерность матрицы очень большая, а нулей в матрице значительно больше ненулевых значений.

----- добавлено через ~11 мин. -----
Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Не понимает, что ноль на диагонали - это автоматическая связь, а не ошибка расчета. Можно ему объяснить, что он неправ?
Нельзя. Матрица становится вырожденной. Нужно или уменьшать размерность такой матрицы исключая строку и стролбец на пересечении нулевого диагонального элемента, либо намеренно изменять этот нулевой элемент на большое число. Сейчас уже точно не помню как это делается в деталях.

----- добавлено через ~13 мин. -----
Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Забавно, что сборка матрицы шла в 10 раз дольше, чем решение системы;
Ну это же си++. Неоптимизированная сборка матрицы и оптимизированное решение системы.

Последний раз редактировалось румата, 25.10.2021 в 09:11.
румата вне форума  
 
Непрочитано 25.10.2021, 10:27
#326
Нубий-IV

Инженер-философ
 
Регистрация: 24.04.2019
Хабаровск
Сообщений: 1,041


Цитата:
Сообщение от nickname2019 Посмотреть сообщение
Нужно в цикле проверять нажатие Esc.
Проблема в том, что при заполнении матрицы и при решении системы циклы не наши, и в них никто ничего не проверяет - ни нули, ни эсекейпы. Либо там есть возможность какие-то свои обработчики на эти события ставить, либо программа будет висеть, пока у пользователя терпение не лопнет. Может, у них это через наследование классов сделано, может - через коллбэки, а может, и вообще невозможно. Пока это еще один открытый вопрос, надо искать ответ.
Цитата:
Сообщение от nickname2019 Посмотреть сообщение
Этот код можно забрать в проект?
Там могут быть ляпы, я его на правильность не проверял, только на скорость.

В первой функции формулы взяты из книги Метод конечных элементов в проектировании транспортных сооружений.А.С. Городецкий, В.И. Зоворицкий, А.И. Лантух-Лященко, А.О. Рассказов. 1981 г..

В остальных функциях формулы получены в Maxima:
В моей библиотеке по МКЭ сорок книг, но полного набора готовых формул я там не нашел. Как вводить шарниры или эксцентриситеты? Как сделать переменный модуль упругости по всем узлам? Какой вектор узловых усилий под неравномерной нагрузкой на части площади КЭ? Какая матрица масс или геометрических параметров у анизотропной плиты? А у гибридного элемента? И т.д. и т.п. Все равно эти формулы придется выводить и тестировать. Вот я заготовку для Максимы и делаю, пока только до балки-стенки дошел.

Думаю, сейчас вообще можно матрицу условными единицами забивать, пока первую версию, работающую от начала и до конца, не получим. По дороге еще придется десяток-другой раз все переделать, потому что никто заранее не предскажет, как сделать правильно. Я бы пока не тратил время на мелочи типа "правильная матрица", пусть там пока везде заглушки стоят.
Цитата:
Сообщение от nickname2019 Посмотреть сообщение
При формировании глобальной матрицы не следует нумеровать степени свободы, на которые наложены связи.
Я имел ввиду связи, забытые пользователем. У Лироскадов есть фишка - обнаруженную в процессе счета изменяемость любезно закреплять автоматической связью, и продолжать расчет. Вот у Старка такой фишки нет - он до первой ошибки считает. И, если таких ошибок полсотни в проекте, придется полсотни раз расчет заново запускать, чтобы их по одному узлу отлавливать.

Eigen при появлении нуля на диагонали просто возвращает в Info код ошибки, а вектор ответа забивает большими числами (кому интересно - файле теста системы размер матрицы уберите в 10, чтобы на экран помещалась, и раскомментируйте строки 36 и 59-62). Ни усилия в схеме посмотреть, ни даже узнать, в каком узле сбой был.
Цитата:
Сообщение от nickname2019 Посмотреть сообщение
лучше использовать полностью-заполненные матрицы, при этом выделять память под локальную матрицу один раз
Это тоже надо тестировать. Может, небольшие матрицы в стеке хранятся, тогда разницы в скорости с C-массивами быть не должно. Например, у меня в одном из тестов вычисление матричное вычисление полинома обгоняло подсчет формулой (хотя и проигрывало схеме Горнера, но формулы матриц вряд ли Горнером сильно оптимизируются).

В книгах по МКЭ много пишут про численное интегрирование. Родной для GMSH решатель GetDP, например, не матрицы элементов хранит, а запрашивает у пользователя число точек Гаусса и типы производных, а остальную арифметику крутит сам - т.е. у него все на матричных произведениях, а не на формулах. Я, когда тест скорости матриц затевал, думал - будет единственный безоговорочный победитель и явные аутсайдеры. А ответ меня удивил сильно. Теперь я никому не верю и все думаю проверять.
Цитата:
Сообщение от румата Посмотреть сообщение
Неоптимизированная сборка матрицы и оптимизированное решение системы.
Обе операции - библиотечные функции. Вот эта троица работает - push_back(), setFromTriplets() и solve(). Все не мои, и никак от меня не зависят. Максимум могу push_back() в восемь разных потоков разогнать, а остальное на совести авторов библиотеки. Тоже открытый вопрос - можно ли полностью параллельное заполнение организовать.
Миниатюры
Нажмите на изображение для увеличения
Название: Maxima.png
Просмотров: 232
Размер:	42.0 Кб
ID:	242046  
Нубий-IV вне форума  
 
Непрочитано 25.10.2021, 11:53
#327
румата


 
Регистрация: 06.04.2015
Сообщений: 1,515


Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Тоже открытый вопрос - можно ли полностью параллельное заполнение организовать.
Из-под .net точно можно, а на ++ я вообще не умею праллелить.
румата вне форума  
 
Автор темы   Непрочитано 25.10.2021, 11:54
#328
nickname2019


 
Регистрация: 18.11.2019
Сообщений: 1,039


Годная книга. На ее основе уже можно линейный процессор написать.
Я раньше пользовался книгой Хечумов Р.А., Кепплер Х., Прокопьев В.И. "Применение метода конечных элементов к расчету строительных конструкций". Более современная книжка "обо всем" в МКЭ, но с опечатками. Найти ее в свободном доступе не могу. Покупать пока не будем, может быть и так справимся.

Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Как вводить шарниры
насколько я понимаю, шарнир - когда убираются линейные степени свободы для данного узла из глобальной системы (так как они закреплены), а угловые степени свободы остаются. Т.е. учет происходит при формировании глобальной матрицы.

Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Как сделать переменный модуль упругости по всем узлам? Какой вектор узловых усилий под неравномерной нагрузкой на части площади КЭ?
Проще всего считать, что во всем КЭ модуль упругости будет постоянный. Также можно усреднить нагрузку для пластин. Для стержней - придется забить формулы, так как балку могут замоделировать одним большим КЭ.

Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Какая матрица масс или геометрических параметров у анизотропной плиты?
Просто матрица D у анизотропной плиты полностью заполнена и получается умножением ортотропной матрицы на геометрическую матрицу из направляющих косинусов. У меня есть материальчик, попробую отсканировать.

Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
А у гибридного элемента? И т.д. и т.п.
Может быть ограничиться оболочечными КЭ (в общем случае анизотропными) и стержневыми?

Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Я имел ввиду связи, забытые пользователем. У Лироскадов есть фишка - обнаруженную в процессе счета изменяемость любезно закреплять автоматической связью, и продолжать расчет. Вот у Старка такой фишки нет - он до первой ошибки считает. И, если таких ошибок полсотни в проекте, придется полсотни раз расчет заново запускать, чтобы их по одному узлу отлавливать.
Да. Надо подумать, можно ли проверить положительную определенность матрицы жесткости до начала длительного расчета. Я пока не знаю как.
nickname2019 вне форума  
 
Непрочитано 25.10.2021, 13:29
1 | #329
RomanM

инженер
 
Регистрация: 25.04.2006
Москва
Сообщений: 1,186


Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Обе операции - библиотечные функции. Вот эта троица работает - push_back(), setFromTriplets() и solve(). Все не мои, и никак от меня не зависят. Максимум могу push_back() в восемь разных потоков разогнать, а остальное на совести авторов библиотеки. Тоже открытый вопрос - можно ли полностью параллельное заполнение организовать.
Не на переаллокацию ли памяти в vector-е всё время сборки матрицы уходит?
Если зарезервировать заранее память (хотя бы приблизительно) перед push_back-ами, должно ускориться.
Дальше возникает вопрос, как бы вообще от ненужного промежуточного vector-а уйти?

...
поэкспериментировал с Вашим кодом, мысли такие.
Во-первых, в релизе соотношение таймингов получается другое - ещё хуже для этапа формирования матрицы. Правда, в абсолютных значениях всё куда быстрее debug.
Далее, думаю, что заполнение SparseMatrix надо анализировать отдельно от заполнения вектора и предыдущих манипуляций.
Там основное время в Вашем коде уходит на арифметику в теле цикла и перевыделение памяти при заполнении вектора. Если это всё убрать (арифметика все равно на практике будет другая), то данный этап ускорится чуть ли не на порядок. Но это просто к слову.
Заполнение SparseMatrix тоже не быстрое. Видимо, основное время уходит на сортировку вставляемых данных.
Можно, кстати, вставлять данные и поштучно непосредственно в SparseMatrix (без промежуточного вектора) через insert,
но, видимо, затраты на пересортировку после каждой вставки могут быть о-очень большие. Так что вариант с передачей сразу всех данных кажется предпочтительным. Хотя, в данной конкретной задаче, если предварительно зарезервировать память в SparseMatrix и выкинуть вектор, затраты времени на эти два этапа у меня сократились на ~40%.

Код:
[Выделить все]
#include <iostream>
#include <ctime>

#include <stdio.h>
#include <tchar.h>
#include <math.h>

#include <Eigen/Dense>
#include <Eigen/Sparse>

int main() {
	int MatrixSize = 716229;
	int BandWidth = 333;

	std::cout << "vector initialization...\n";
	unsigned start_time = std::clock();
	std::vector<Eigen::Triplet<double>> T;
	T.reserve(MatrixSize*BandWidth);
	for (int i = 0; i < MatrixSize; i++) {
		int bandMinIndex = i;
		int bandMaxIndex = std::min(i + BandWidth, MatrixSize);
		for (int j = bandMinIndex; j < bandMaxIndex; j++) {
			//double distance = i - j;
			//double max_distance = BandWidth;
			//double relativeDistance = distance / max_distance;
			//T.push_back(Eigen::Triplet<double>(i, j, 1.0 - std::pow(relativeDistance, 2)));
			T.push_back(Eigen::Triplet<double>(i, j, 1.));
			//K.insert(i, j) = 1.0 - std::pow(relativeDistance, 2);
		}
	}
	std::cout << (std::clock() - start_time) / 1000.0 << " seconds\n" << "...done.\n\n";
	
	std::cout << "SparseMatrix initialization...\n";
	start_time = std::clock();
	Eigen::SparseMatrix<double> K(MatrixSize, MatrixSize);
	//K.reserve(Eigen::VectorXi::Constant(MatrixSize,BandWidth));
	K.setFromTriplets(T.begin(), T.end());
	std::cout << (std::clock() - start_time) / 1000.0 << " seconds\n" << "...done.\n\n";

	std::cout << "Forces vector calulation...\n";
	start_time = std::clock();
	Eigen::VectorXd F(MatrixSize);

	for (int i = 0; i < MatrixSize; i++)
		F[i] = (double)i / MatrixSize;

	unsigned end_time = std::clock();
	std::cout << (end_time - start_time) / 1000.0 << " seconds\n";
	std::cout << "...done.\n\n";

	std::cout << "Equation system solving...\n";
	start_time = std::clock();

	Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>> solver(K);
	Eigen::VectorXd x = solver.solve(F);

	end_time = std::clock();
	std::cout << (end_time - start_time) / 1000.0 << " seconds\n";
	std::cout << "...done.\n";

	system("pause");

	//std::cout << "Solver info: " << solver.info() << std::endl;
	//std::cout << "\nStiffness:\n" << K << std::endl;
	//std::cout << "\nForces:\n" << F << std::endl;
	//std::cout << "\nDisplacements:\n" << x << std::endl;

	return 0;
}
Миниатюры
Нажмите на изображение для увеличения
Название: 2021-10-26_00-09-04.png
Просмотров: 19
Размер:	5.1 Кб
ID:	242078  

Последний раз редактировалось RomanM, 26.10.2021 в 00:44.
RomanM вне форума  
 
Непрочитано 25.10.2021, 14:07
1 | #330
Бахил

?
 
Регистрация: 17.06.2014
Царицын
Сообщений: 10,335


Цитата:
Сообщение от nickname2019 Посмотреть сообщение
можно ли проверить положительную определенность матрицы жесткости до начала длительного расчета
Можно проверить на наличие нулей на главной диагонали. И, либо выдать предупреждение, либо вставить любое число, лучше 1е10.

----- добавлено через ~2 мин. -----
Есть Холецкий на Паскале. Если интересно, могу скинуть. Правда только однопроцессорный.
__________________
В конструктивных дискуссиях каждый участник укрепляется в своих заблуждениях.
Бахил вне форума  
 
Автор темы   Непрочитано 25.10.2021, 14:26
#331
nickname2019


 
Регистрация: 18.11.2019
Сообщений: 1,039


Цитата:
Сообщение от Бахил Посмотреть сообщение
Есть Холецкий на Паскале. Если интересно, могу скинуть. Правда только однопроцессорный.
Холецкий должен быть в eigene.

Последний раз редактировалось nickname2019, 25.10.2021 в 15:00.
nickname2019 вне форума  
 
Непрочитано 26.10.2021, 05:33
#332
Нубий-IV

Инженер-философ
 
Регистрация: 24.04.2019
Хабаровск
Сообщений: 1,041


Цитата:
Сообщение от nickname2019 Посмотреть сообщение
Проще всего считать, что во всем КЭ модуль упругости будет постоянный.
Раз заявлена нелинейность для железобетона - скорее всего придется интерполяцию жесткостей делать, вдобавок к анизотропии, иначе придется мельчить сетку. Формально это не сложно - достаточно в файле максимы добавить строку типа E=E1+(E2-E1)u + (E3-E1)v, и посмотреть ответы.
Цитата:
Сообщение от nickname2019 Посмотреть сообщение
Может быть ограничиться оболочечными КЭ
Даже оболочек - и то целый выводок в любой книге. Такие полиномы - сякие полиномы - а вот кому косинусы. Помнится, в старой Лире наклонные оболочки были геометрически изменяемыми, а в новых версиях - уже нет; явно там в элементах что-то подшаманили (или в автоматических связях, кто их там разберет, закрыто-коммерческих). Даже из Зенкевича полно ссылок на конкретные элементы, в виде научных статей, за которые денег хотят.

Мне вот Старковские элементы нравятся - и повышенный порядок, и узловые усилия, и плоскости осреднения. Старк благодаря ним даже с сеткой не заморачивается - все прямоугольным шаблоном делает, и только на стыках типовую разбивку примыкающих элементов применяет. Но формул на них нет в свободном доступе, тут только подбирать полиномы и тестировать.

В любом случае полный набор формул "для всего" вряд ли хоть в одной книге найдем; что-то придется выводить. Во вложении - пример матрицы в Максиме, и pdf с описанием процесса.
Цитата:
Сообщение от румата Посмотреть сообщение
Из-под .net точно можно
В плюсах многопоточность как часть языка сделали стандартом только в последних редакциях, из VS2008 она недоступна. Но через WinAPI можно это сделать. Вопрос в другом - поддерживает ли сама библиотека многопоточность. Решатели-то вроде все ядра грузят. А вот доступ к матрице на добавление элементов - не факт, что можно в несколько потоков пустить: если про это не подумали разработчки, есть шанс просто запортить внутренние указатели.
Цитата:
Сообщение от Бахил Посмотреть сообщение
Можно проверить на наличие нулей на главной диагонали.
Не поможет. Пока решено нулевые элементы в глобальную матрицу не записывать, так что случая "до расчета все нули в строке" у нас не предполагается. А если пользователь просто связь забыл - то до расчета на диагонали нуля нет, он всплывет только в процессе решения, когда по всем узлам недозакрепленной части конструкции расчет пройдет - только на последнем из них этот ноль и выскочит.

По-хорошему, у решателя должен быть способ установить обработчик такого события. Некая функция, которая вызывается после обнаружения нуля, выполняет нужные действия (ставит связь, пишет протокол), возвращает решателю свой ОК, и тот продолжает счет. Сходу я такого не нагуглил.
Вложения
Тип файла: zip Балка-стенка.zip (712.9 Кб, 15 просмотров)
Нубий-IV вне форума  
 
Непрочитано 26.10.2021, 07:44
#333
Бахил

?
 
Регистрация: 17.06.2014
Царицын
Сообщений: 10,335


Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Вложения
А это что за опус? Offtop: Предупреждать же надо.
__________________
В конструктивных дискуссиях каждый участник укрепляется в своих заблуждениях.
Бахил вне форума  
 
Автор темы   Непрочитано 26.10.2021, 08:20
#334
nickname2019


 
Регистрация: 18.11.2019
Сообщений: 1,039


Вот вывод МЖ прямоуголных КЭ плоской задачи и физ. нелинейного пространственного КЭ - параллелепипеда.

Это для справки по учету нелинейности в МЖ.
Вложения
Тип файла: pdf Вывод МЖ прямоугольных КЭ сплошного тела.pdf (1.43 Мб, 24 просмотров)
nickname2019 вне форума  
 
Непрочитано 26.10.2021, 09:20
#335
румата


 
Регистрация: 06.04.2015
Сообщений: 1,515


Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Не поможет. Пока решено нулевые элементы в глобальную матрицу не записывать, так что случая "до расчета все нули в строке" у нас не предполагается.
Если нулевые элементы в матрицу не записывать, то получится разреженная матрица. Не вижу проблемы в проверке на диагональный ноль-элемент в разреженной матрице. Кроме того, возможно, еще прийдется формировать матрицу индексов, в ней точно нулевые элементы на диагонали будут видны.
Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
По-хорошему, у решателя должен быть способ установить обработчик такого события.
Ну можно попробовать и через функцию обратного вызова решить. Возможно так будет отыскание нулей на диагонали более рациональным.
Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Даже оболочек - и то целый выводок в любой книге. Такие полиномы - сякие полиномы - а вот кому косинусы.
Без этого обойтись не получится. Все эти балки-стенки и исключительно прямоугольные оболочки никому не нужны. Нужны произвольные четырехугольные и треугольные КЭ оболочки. Потом четырехугольные и треугольные КЭ толстой оболочки или плиты.
румата вне форума  
 
Автор темы   Непрочитано 26.10.2021, 09:49
#336
nickname2019


 
Регистрация: 18.11.2019
Сообщений: 1,039


Цитата:
Сообщение от румата Посмотреть сообщение
Если нулевые элементы в матрицу не записывать, то получится разреженная матрица. Не вижу проблемы в проверке на диагональный ноль-элемент в разреженной матрице. Кроме того, возможно, еще прийдется формировать матрицу индексов, в ней точно нулевые элементы на диагонали будут видны.
На самом деле, "истинно" нулевых элементов на диагонали матрицы не будет. Будут "околонулевые" элементы. Возможно, имеет смысл на этапе формирования матрицы выичислять среднее значение диагонального элемента. Если есть элементы, которые на порядки меньше среднего значения -подсвечивать эти узлы в схеме или элементы, примыкающие к этим узлам (или сразу задавать большое число на диагонали матрицы в этом месте, как предлагал Бахил).

Цитата:
Сообщение от румата Посмотреть сообщение
Все эти балки-стенки и исключительно прямоугольные оболочки никому не нужны. Нужны произвольные четырехугольные и треугольные КЭ оболочки. Потом четырехугольные и треугольные КЭ толстой оболочки или плиты.
Есть вариант вывести матрицу жесткости пространственного КЭ (см. вложение пост #334 рис. 3.5, но не прямоугольный, а произвольный), а из нее получить КЭ оболочки (единичное смещение давать не одному узлу, а сразу двум, потом получать реакции, вычислять суммарную реакцию - переводить в M,N,Q).
Может быть, это будет сложнее, а может проще, не могу пока понять.
nickname2019 вне форума  
 
Непрочитано 26.10.2021, 10:18
#337
Нубий-IV

Инженер-философ
 
Регистрация: 24.04.2019
Хабаровск
Сообщений: 1,041


Цитата:
Сообщение от румата Посмотреть сообщение
Не вижу проблемы в проверке на диагональный ноль-элемент в разреженной матрице.
А она есть .
Элемент фермы в плоской задаче со степенями свободы x,y.
  1. Заносим элемент в матрицу жесткости. Поскольку по оси y реакций элемент не передает, вторая и последняя строки (суммы сил на вертикаль в 1 и 2 узле) - полностью нулевые.
  2. Удаляем строки с нулями на диагонали (так делают Лира и Скад, но не делает Старк). Это, фактически, значит "Ставим автоматические связи".
  3. Начинаем расчет. Исключаем из второй строки первый элемент: добавляем в нее первую строку.
  4. Опа! Нарисовался! Ноль на диагонали. До начала расчета его там не было. Тут Лира и Скад ставят на диагональ ненулевое значение ("Автоматическая связь в процессе расчета"), а Старк падает в обморок. У Eigen к этому месту нет доступа - solve() вернет ошибку, и абсолютно левые перемещения.

Цитата:
Сообщение от румата Посмотреть сообщение
Все эти балки-стенки и исключительно прямоугольные оболочки никому не нужны.
Как минимум они нужны авторам книг, чтобы показать, какие простые формулы получаются в МКЭ:
Насколько я понимаю, "оболочка" = "балка-стенка" + "плита". Плита у меня следующая в списке. А дальше можно будет анизотропию добавить; и, глядя на эти формулы, хочется проверить, не будет ли там численное интегрирование быстрее.
Миниатюры
Нажмите на изображение для увеличения
Название: Изменяемость.png
Просмотров: 140
Размер:	20.7 Кб
ID:	242090  Нажмите на изображение для увеличения
Название: Элемент.png
Просмотров: 138
Размер:	122.0 Кб
ID:	242093  
Нубий-IV вне форума  
 
Автор темы   Непрочитано 26.10.2021, 10:34
#338
nickname2019


 
Регистрация: 18.11.2019
Сообщений: 1,039


Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Плита у меня следующая в списке. А дальше можно будет анизотропию добавить; и, глядя на эти формулы, хочется проверить, не будет ли там численное интегрирование быстрее.
Я выводил матрицу прямоугольного пространственного анизотропного КЭ 24х24 элемента. Получилась простыня.

Потом понял, что это бессмысленно (когда решил проверить полученную матрицу численным интегрированием). Быстрее и проще проинтегрировать численно с использованием весовых коэффициентов. Объем вычислений тот же.
Если степень полинома известна, то можно так подобрать точки для вычисления интеграла, что результаты интегрирования будут точными.
nickname2019 вне форума  
 
Непрочитано 26.10.2021, 10:58
#339
румата


 
Регистрация: 06.04.2015
Сообщений: 1,515


Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
и, глядя на эти формулы, хочется проверить, не будет ли там численное интегрирование быстрее.
Однозначно будет быстрее. На порядки. Все матричные вычисления нужно делать численно на числах с двойной точностью. А вся эта символьная борода может быть полезной только для наглядности.
Цитата:
Сообщение от Нубий-IV Посмотреть сообщение
Элемент фермы в плоской задаче со степенями свободы x,y.
1)Вычисляется матрица жесткости элемента Kr в локальной системе координат.
2)Вычисляется матрицы преобразования матрицы жескости элемента из локальной в глобальную систему координат Tr
3)Учитываеются элементные шарниры путем вычисления матрицы L и последующего преобразования матрицы Кr = L^T*Kr*L
4)Вычисляется матрица жесткости элемента в глобальной системе координат К'r = Tr^T*Kr*Tr
5)Формируется общая матрица жесткости системы К из матриц элементов К'r по матрице индексов или иным способом суммирования.
6)Выполняется проверка матрицы К на наличие нулевых элементов на диагонали перед решением СЛАУ или обращением этой матрицы.
Где эта проблема, которая есть?
румата вне форума  
 
Автор темы   Непрочитано 26.10.2021, 11:12
#340
nickname2019


 
Регистрация: 18.11.2019
Сообщений: 1,039


Цитата:
Сообщение от румата Посмотреть сообщение
Однозначно будет быстрее. На порядки. Все матричные вычисления нужно делать численно на числах с двойной точностью. А вся эта символьная борода может быть полезной только для наглядности.

1)Вычисляется матрица жесткости элемента Kr в локальной системе координат.
2)Вычисляется матрицы преобразования матрицы жескости элемента из локальной в глобальную систему координат Tr
3)Учитываеются элементные шарниры путем вычисления матрицы L и последующего преобразования матрицы Кr = L^T*Kr*L
4)Вычисляется матрица жесткости элемента в глобальной системе координат К'r = Tr^T*Kr*Tr
5)Формируется общая матрица жесткости системы К из матриц элементов К'r по матрице индексов или иным способом суммирования.
6)Выполняется проверка матрицы К на наличие нулевых элементов на диагонали перед решением СЛАУ или обращением этой матрицы.
Где эта проблема, которая есть?
Если юзер забыл приложить связи на систему, то в процессе решения может быть ошибка (если матрицу решать гауссом, то на диагонали в процессе решения появиться ноль).
nickname2019 вне форума  
Ответ
Вернуться   Форум DWG.RU > Программное обеспечение > Программирование > Как создать расчетное программное обеспечение с открытым исходным кодом (конструктивные решения)

Размещение рекламы
Опции темы Поиск в этой теме
Поиск в этой теме:

Расширенный поиск


Похожие темы
Тема Автор Раздел Ответов Последнее сообщение
СП 335.1325800.2017 «Крупнопанельные конструктивные системы. Правила проектирования» (Обсуждение) Armin Прочее. Архитектура и строительство 37 07.11.2018 06:55
Фирменные решения по пропуску коммуникаций через стены подвала Regby Конструкции зданий и сооружений 2 07.04.2010 20:43
устройство и возможные конструктивные решения вентфасада из кирпича Ivansobaka Каменные и армокаменные конструкции 1 16.12.2009 06:38
Конструктивные решения по перемычкам в многослойных кирпичных стенах! Westroy Архитектура 16 30.11.2009 13:57
Конструктивные решения монтажных соединений многоэтажных зданий на высокопрочных болтах VoRoNoFF Конструкции зданий и сооружений 1 04.04.2009 00:41