Comments 46
С барьерами памяти есть одна небольшая подстава: на каждой архитектуре они работают по-разному. Так, для x86 справедливо следующее:
Loads and stores are not reordered with locked instructions
То есть любая атомарная операция на x86, кроме load/store, которые представляют собой обычный MOV — это неявный барьер. Отсюда получается, что в многопоточном коде под x86 барьеры практически никогда не используются, а компилятор C++ генерирует одинаковый код вне зависимости от используемой модели. Бывают, конечно, редкие исключения: запись в память, минуя кэш, или же синхронизация чисто на load/store без блокировки шины. Но это всё экзотика.
К чему я клоню: под x86 можно годами писать параллельный код, вообще не зная о существовании барьеров. А вот с GPU, увы, уже так не получится.
Выяснил, что на ARM'ах, с которыми приходится работать уже на основной работе, тоже есть опасность подобных ситуаций. К счастью, стандарт C++ автоматически защищает несведущих людей, используя по-умолчанию медленный, но безопасный std::memory_order_seq_cst.
Два доклада Фёдора Пикуса многое прояснили в этой области. Поэтому я решил ими поделиться.
К счастью, стандарт C++ автоматически защищает несведущих людей, используя по-умолчанию медленный, но безопасный std::memory_order_seq_cst.
Он защищает только при использовании типов std::atomic<>
. Берём любую самодеятельность: да хоть atomicExch
или же interlocked API над интовыми указателями — и всё ломается.
Целью же является проверка правдоподобности модели подъёмной силы Ньютона. А именно, будут ли наблюдаться какие-либо интересные эффекты, если взять большое количество шариков и профиль крыла.
на основе шариков и плавают, примерно на тех же принципах как Ваше крыло летит )
Красота картинки важна.
Способность нарисовать картинку, похожую на реальность, она похвальна и осмысленна — в тех же компьютерных играх, к примеру. Но физика, вроде бы, все же претендует на описание настоящей реальности…
Температурная карта давленияЭто финиш… Она температурная или давления? Какая такая карта, у вас просто картина распределения в пространстве. Нормально это называть можно так: «Контуры давления», «Распределение поля давления», «Распеределение давления», «Поле давления». Легенды где?
Один из многих способов численно решить систему уравнений — это применить Метод Конечных Элементов в связке с Адаптивным Сеточным Методом в случае подхода Эйлера.Не вводите читателей в заблуждение, основным, самым используемым способо решения является метод конечных объемов. SpaceX конкретно не раскрыли что они использовали, но судя по всему это был метод конечеых объемов или конечно разностный метод.
Эта система определяет равнодействующую всех сил (давления, вязкости и гравитации) для всех направлений. Дополнительно может быть задано второе векторное уравнение, обеспечивающее условие несжимаемости, если нужно описать жидкости
Только одно уравнение из уравненИЙ Навье-Стокса это уравнение баланса импульса. «Дополнительно», дополнительно блин, уравнение неразрывности с каких пор у нас стало дополнительным?! При чем здесь описание несжимаемости?! Как вы замыкаете эти два уравнения? Где уравнение состояния жидкости?
Для иных сред, чаще всего реализуют подход Лагранжа, через метод сглаженных частиц (Smoothed-particle hydrodynamics, SPH).Да, чаще всего для рисования красивых картинок, не более того, у него огромные проблемы с постановкой корректных ГУ и с сохранением массы. Все эти качественные красивые эффекты (брызги пена и тд) туфта, когда они появляются не четко тогда, там и таким образом как в реальном процессе.
Для приемлимых результатов с Лагранжевым подходом используют метод FPM — метод конечных точек.
И самое интересное, почему в явно сверхзвуковом режиме течения вокруг крылового профиля у вас есть возмущения выше по потоку в расчете? Какие скорости течения вы использовали, как замыкали уравнения навье стокса? Какой жидкости соответствует ваша модель?
Пока на физику это не похоже, на красивые мультики да.
Самый главный вопрос, где подъемная сила?!
Несоответствие заголовка и содержания. Я зашел почитать про расчет подъемной силы!
Назовите статью тогда «Рисование красивых картинок навеянных физикой методом частиц на CUDA»
«Интересные эффекты», визуальные эффекты?
А какие вы параметры брали? Соотношение скорости звука и скорости крыла и размера крыла точно реалистичные? Потому что сверх звуковой полёт это уже другое.
Проблема в том, что устойчивое название такой палитры — "тепловая карта", а не "температурная карта".
Проблема только в том, что за общепринятый стандарт начали выдавать какой-то бред, вылезший из каких-то людей видевших то ли изображения с тепловизора то ли расчеты теплообмена и решивших, что это всё про теплоту, даже не температуру.
А про «радугу» впервые слышу, в ваших же перечисленных пакетах, также используется понятие heatmap: www.ansys.com/blog/heat-maps-determine-build-orientations-additive-manufacturing
«Температурная карта» — это вроде устойчивое название такой палитры…«Палитры», понимаете палитры, палитра называется «радуга».
Высосаный из пальца источник, где heat map применяется для описания отклика от двух параметрво угла наклона. В практике ANSYS никогда не используется heat map для описания пространственных распределений параметров, нате вам ссылку тоже: www.ansys.com/products/fluids/ansys-cfx просто посмотрите как называется распределение давления на картинке.
Я против размывания вполне конкретных устоявшихся терминов обеспечивающих максимальную корректность. Я вот работаю в сфере гидрогазодинамики и термодинамики, и если мои коллеги по цеху начнут писать «heat flux heat map», то я очень расстроюсь, каша в словах — каша в головах. Автор, очевидно, не специалист в CFD, и не знаком с физикой, а вы его защищаете, зачем?
Однако, у этой структуры есть следующий ряд фундаментальных ограничений:
1. Множество значений координат ограничено размером самой сетки.
нет, не ограничено, для бесконечного пространства координата приводится к ячейке как остаток от деления на ширину/высоту сетки, играя роль хеша, а при обработке столкновений дополнительно проверяется реальное расстояние между парой частиц-кандидатов
насчет остальных аргументов — большинство также вызывают скепсис:
интересно сравнить конкретные реализации по производительности
вы пишите про 8-12 FPS на 1080 Ti на 2млн частиц, у меня получатся считать на ~50-60 FPS это количество на таком железе, можете попробовать запустить демо-версию, там есть практически 1:1 такая же простая симуляция-коллайдер:
store.steampowered.com/app/1196080/Space_Simulation_Toolkit
работает на основе inclusive_scan, без fast_radix_sort, с обычной сеткой, это самый быстрый способ из всех, что я исследовал:
on-demand.gputechconf.com/gtc/2014/presentations/S4117-fast-fixed-radius-nearest-neighbor-gpu.pdf
применение вообще любой сортировки на каждом кадре симуляции делает построение сетки самой дорогой операцией (примерно x10 по отношению ко всему остальному коду симуляции)
именно поэтому fluids_v3 алгоритм такой крутой, сетку можно строить вообще без единой сортировки, за счет prefix_sum (inclusive_scan)
вот 400K частиц на 780 Ti, даже без VBO_INTEROP оптимизации рендера
вы пишите про 8-12 FPS на 1080 Ti на 2млн частиц, у меня получатся считать на ~50-60 FPS это количество на таком железе
Тут довольно тонкий момент. Процитирую один фрагмент:
Модули слабо связаны друг с другом и напрямую не зависят. Существует возможность проводить симуляцию в оффлайн режиме и рендерить результаты позднее. В CUDA-версии демонстрации рендеринг, что на видео в начале статьи, происходит каждые 16 шагов.
То есть эти 12 FPS нужно расценивать как 192 шага в секунду. Причина по которой я рисую не каждый шаг, а через 16 шагов в том, что передача данных между контекстами CUDA и OpenGL занимает существенное время. Я тем самым уплотняю вычисления и минимизирую задержки на передачу данных через PCI-E шину.
Мне известно о наличии OpenGL Interoperability в CUDA API, но конкретно в этом проекте я хотел получить возможность рендерить на другом видеоадаптере. Поэтому модуль, занимающийся расчётами, хоть и использует видеоадаптер для вычислений, никак не использует его для отрисовки.
Если вдруг интересно, предлагаю ознакомиться с вот этим файлом:
github.com/tony-space/WingSimulator/blob/master/Simulation/libCudaSimulation/CSimulationCuda.cu
А именно, вот этот фрагмент:
float CSimulationCuda::Update(float dt)
{
for (int i = 0; i < 16; ++i)
{
m_odeSolver->NextState(ComputeMinDeltaTime(dt), m_curOdeState, m_nextOdeState);
m_nextOdeState.swap(m_curOdeState);
}
ColorParticles2();
return dt;
}
Основной недостаток этой модели — это отсутствие взаимодействия частиц газа друг с другом
Это потому что межатомное взаимодействие проявляется на расстояниях до 6-8 ангстрем, а при комнатной температуре и атмосферном давлении среднее расстояние между частицами газа порядка 30 ангстрем.
У меня некоторый диссонанс остался, вначале вы пишите о симуляции из отдельных атомов, а потом говорите о конечных элементах и прочих вещах, где единицей является какая-то более обобщенная структура. Но за статью спасибо, узнал кое-что интересное, например о match_any, может быть удастся применить с пользой.
А как выглядит этим методом дозвуковой режим?
Фраза
Основной недостаток этой модели — это отсутствие взаимодействия частиц газа друг с другом. Из-за этого при нормальных условиях она даёт некорректные результаты, хотя всё ещё может применяться для экстремальных условий, где взаимодействием можно пренебречь.
прямо пугает. На мой взгляд, модель идеального газа — это и есть нормальные условия, ничего экстремального.
Тут на первой же картинке очевидно, что дело идет о сверхзвуковом обтекании — ясно видна ударная волна. Судя по углу ее конуса (на глаз), там почти 3M. Так что частицам приходится взаимодействовать.
А как выглядит этим методом дозвуковой режим?
Не стоит воспринимать модель подъёмной силы Ньютона (Newtonian Lift) как адекватную модель для точных газодинамических расчётов. На сайте NASA, ссылка в аннотации перед катом, подробно описано почему она некорректна. Это всего лишь грубая модель взаимодействия газа с профилем крыла.
На мой взгляд, модель идеального газа — это и есть нормальные условия, ничего экстремального.
Некоторые модели идеального газа допускают эластичные столкновения между частицами. Главное чтобы не было дальнодействующих сил.
en.wikipedia.org/wiki/Ideal_gas#Microscopic_model
Сам думаю добраться до лагранжевых методов, но пока всё руки не доходят.
Позволю пару придирок.
«Метод Конечных Элементов» — вряд ли заслуживает больших букв — и даже Википедия на это намекает. Большие буквы в научно-популярной статье сразу пугают.
По поводу тепловых/температурных карт — что-то ни разу не слышал такого жаргонного употребления (а на конференциях по гидродинамике я регулярно бываю). Все говорят — «поле скорости» или «распределение давления». Возможно это некий внутренний жаргон? Но может я просто из другой песочницы — век живи — век учись.
По поводу уравнений Навье-Стокса. Эта система определяет равнодействующую всех сил (давления, вязкости и гравитации) для всех направлений.
Вообще-то это уравнение сохранения импульса для жидкого объема. Гравитация там слегка лишняя — в том смысле, что там сумма сил куда может входить всё что угодно, например электромагнитные силы.
По поводу МКЭ. Еще в 2007 году, Брайен Сполдинг призывал всех отказаться от МКЭ в пользу МКО. Что в общем-то и происходит — пакеты которые популярны используют МКО. Научная среда — другая история, там у каждого свой код.
По поводу подъемной силы:
«На тепловой карте (heatmap) давления ниже отчётливо видны области высокого и низкого давления.» Ну как бы не очень. Видны цветные области — что они значат — «догадайся мол сама». Но не все читатели регулярно видят такие картинки — желательно пояснить.
Как у вас взаимодействуют частицы? Что происходит при столкновении? Твердые сферы? Возможно я проглядел этот момент.
За урок по Куде — отдельное спасибо.
Из Лагранжевых методов советую finite point method.
Кстати, а какая сейчас ситуация с «гибридными» методами вроде метода крупных частиц Белоцерковского и Давыдова? Используют сейчас или не прижились?
Развитие методов с Лагранжевым подходом идет в задачах со сложной реологией жидкостей (и не совсем жидкостей), где вязкости в одной области могут отличаться на многие порядки, стекло, пластики, паста.
Гибрид Лагранжева и Эйлерова подхода применяется для многофазных течений, для гомогенных сред нет.
Ну, меня не столько метод крупных частиц интересовал, сколько сам гибридный подход. Просто я со времён университета этим не занимаюсь, и он первым на ум пришёл. Да и даже тогда мы его только в учебных целях программировали, а для реальных расчётов, я помню, брали метод Годунова и его модификации вроде PPM. Но то были астрофизические расчёты, там не было каких-то сложных граничных условий (зато было много других проблем).
«Метод Конечных Элементов» — вряд ли заслуживает больших букв — и даже Википедия на это намекает. Большие буквы в научно-популярной статье сразу пугают.
Замечание дельное, исправил. Спасибо.
По поводу тепловых/температурных карт — что-то ни разу не слышал такого жаргонного употребления (а на конференциях по гидродинамике я регулярно бываю). Все говорят — «поле скорости» или «распределение давления». Возможно это некий внутренний жаргон? Но может я просто из другой песочницы — век живи — век учись.
Я знаком с термином прежде всего из англоязычных источников (Heat map). Не исключаю, что в узком профессиональном сообществе приняты другие названия.
По поводу уравнений Навье-Стокса. Эта система определяет равнодействующую всех сил (давления, вязкости и гравитации) для всех направлений.
Вообще-то это уравнение сохранения импульса для жидкого объема. Гравитация там слегка лишняя — в том смысле, что там сумма сил куда может входить всё что угодно, например электромагнитные силы.
Опять же, знаком с уравнениями из англоязычных источников, где они записывается в виде:
. Допускаю, что в русскоязычной литературе принята другая форма. Из того, что я в своё время понял, это просто сумма сил, действующих на элементарный объём dV массой dm. Собственно поэтому у нас и фигурирует плотность как отношение этих двух бесконечно малых величин. Если обе части мы домножим на dV, мы получим силу в Ньютонах, правда бесконечно малую. Поправьте, если я что-то не правильно понял.
По поводу МКЭ. Еще в 2007 году, Брайен Сполдинг призывал всех отказаться от МКЭ в пользу МКО. Что в общем-то и происходит — пакеты которые популярны используют МКО. Научная среда — другая история, там у каждого свой код.
Замечание разумное. В сеточные методы не углублялся, так как больше интересуюсь бессеточными, поэтому вполне могу обладать устаревшей информацией. Просто интересно, вырождается ли в двумерном случае метод конечных объёмов в МКЭ?
«На тепловой карте (heatmap) давления ниже отчётливо видны области высокого и низкого давления.» Ну как бы не очень. Видны цветные области — что они значат — «догадайся мол сама». Но не все читатели регулярно видят такие картинки — желательно пояснить.
Добавил текстовую расшифровку в указанное предложение. Цветовую легенду всё-же не стал прикладывать, потому что в достаточно грубой модели Newtonian Lift значения давления и ускорений являются крайне не точными. Труъ-физикам они покажутся чересчур неадекватными.
Как у вас взаимодействуют частицы? Что происходит при столкновении? Твердые сферы? Возможно я проглядел этот момент.
Частицы ведут себя как неупругие шарики, примерно так же, как и в DEM. В частности демпфирование пришлось добавить для улучшения численной стабильности. Можете посмотреть как используется функция SpringDamper: github.com/tony-space/WingSimulator/blob/master/Simulation/libCudaSimulation/CDerivativeSolver.cu
Кажется я понял в чем дело.
«Тепловая карта» — относительно новый термин.(С) Википедия.
В исследованиях, гидродинамика и теплообмен часто идут в одном флаконе. И на соседних картинках может быть и распределение температуры и скажем распределение давления. И если называть это одним словом то… председатель секции (или рецензент) будет недоволен.
Насчет Навье-Стокса — суть-то не в силах, а в скорости которую можно извлечь из ускорения. В каком нибудь Навье-Стоксе редуцированном до уравнения Эйлера вообще с силами грустно. Ну я думаю, вы всё правильно понимаете — другой вопрос, что написать популярную статью подобрав точные фразы для всего вовсе не просто.
Так, а что такое "(скалярная сумма модулей сил)"?. Давление это сила делённая на площадь. Или сила — градиент давления. У вас должна быть средняя сила за какой-то период осреднения отнесенная к элементарной площадке — что касается распределения давления по крыловому профилю.
Если про диссипацию и потенциал взаимодействия напишете подробнее — будет интересно. В той программе я вряд ли разберусь. Собственно, вся картина и определяется же потенциалом взаимодействия и количеством частиц.
Опять же, знаком с уравнениями из англоязычных источников, где они записывается в виде:
image. Допускаю, что в русскоязычной литературе принята другая форма.
Это шутка такая что-ли? В любой, абсолютно любой (русско-, англо-, немецко-, китайско- итд язычной) литературе минимальный набор уравнений Навье-Стокса представлен двумя дифференциальными уравнениями сохранения импульса и сохранения массы или неразрывности (без разложения по направлениям)! Вы тут пишите одно. У меня вопрос, если вы не способны понять английский текст даже в той же википедии, может вам по русски всё же читать с начала?
Ошибка «Run-Time Check Failure #3 — The variable 'result' is being used without being initialized.»
CUDA версии 11.0, VS 2019. Сэмплы CUDA собираются и работают «искаропки».
libCudaSimulation.dll!thrust::detail::allocator_system<thrust::device_allocator>::get(thrust::device_allocator & a) Line 360 C++
libCudaSimulation.dll!thrust::detail::allocator_traits_detail::copy_construct_range<thrust::cuda_cub::tag,thrust::device_allocator,thrust::detail::normal_iterator<thrust::device_ptr>,thrust::device_ptr>(thrust::execution_policy<thrust::cuda_cub::tag> & from_system, thrust::device_allocator & a, thrust::detail::normal_iterator<thrust::device_ptr> first, thrust::detail::normal_iterator<thrust::device_ptr> last, thrust::device_ptr result) Line 218 C++
libCudaSimulation.dll!thrust::detail::copy_construct_range<thrust::cuda_cub::tag,thrust::device_allocator,thrust::detail::normal_iterator<thrust::device_ptr>,thrust::device_ptr>(thrust::execution_policy<thrust::cuda_cub::tag> & from_system, thrust::device_allocator & a, thrust::detail::normal_iterator<thrust::device_ptr> first, thrust::detail::normal_iterator<thrust::device_ptr> last, thrust::device_ptr result) Line 291 C++
libCudaSimulation.dll!thrust::detail::contiguous_storage<float2,thrust::device_allocator>::uninitialized_copy<thrust::detail::normal_iterator<thrust::device_ptr>>(thrust::detail::normal_iterator<thrust::device_ptr> first, thrust::detail::normal_iterator<thrust::device_ptr> last, thrust::detail::normal_iterator<thrust::device_ptr> result) Line 284 C++
libCudaSimulation.dll!thrust::detail::vector_base<float2,thrust::device_allocator>::append(unsigned __int64 n) Line 877 C++
libCudaSimulation.dll!thrust::detail::vector_base<float2,thrust::device_allocator>::resize(unsigned __int64 new_size) Line 327 C++
libCudaSimulation.dll!wing2d::simulation::cuda::CSimulationCuda::CopyToGPU() Line 179 C++
libCudaSimulation.dll!wing2d::simulation::cuda::CSimulationCuda::ResetState(const wing2d::simulation::SimulationState & state) Line 113 C++
WingSimulator.exe!SetupState(wing2d::simulation::ISimulation * simulation, std::vector<std::tuple<float,float>,std::allocator<std::tuple<float,float>>> && airfoil) Line 117 C++
WingSimulator.exe!main(int argc, char * * argv) Line 132 C++
Библиотека thrust в очередной раз расстраивает.
P.S. можете зайти в раздел релизов на гитхабе и попробовать запустить бинарник через заранее заготовленные BAT файлы.
Проверить текущую версию можно в «Панель управления NVIDIA» -> «Справка» -> «Информация о системе»:
P.S. А драйвер у меня референсный, версии 27.21.14.5182 от 7/19/2020.
P.S. Если хотите, можем перейти на гитхаб для обсуждения, тут не очень подходящее место.
Насколько я понял статью про подъёмную силу крыла, то там присутствуют такие эффекты:
1. Отрицательного давления не бывает, поэтому верхняя кромка крыла не может тянуть вверх, а только не препятствовать этому, соотвественно, подъёмной силой обладает только нижняя плоскость и как следствие — угол атаки, это решающий фактор.
2. На нижнюю кромку давят частицы не один раз, как в Ньютоновской модели, а переотражаются многократно между крылом и набегающим потоком и фактически, почти все потоки, кроме кромочного проходят мимо, сталкиваясь с подушкой, а не крылом, как было бы в Ньютоновской модели. Но куммулятивный эффект равнозначен движению в более плотной среде, это очень сильно увеличивает эффект подъёмной силы.
3. Если происходит срыв потока, то через переднюю кромку выдавливается против движения подушка в верхнюю часть и переотражения нарушаются, что приводит к резкому падению подъёмной силы. Особой формой выгнутости крыла, можно заполнить область пониженного давления и отсрочить тем самым предавливание потока, при незначительном падении подъёмной силы (если форма верхней кромки будет совпадать по форме с областью низкого давления, то это давление не особо возрастёт), поэтому добиться больших углов атаки без данного эффекта. Поэтому же самолёты легко летают вверх ногами, просто больше рискуют свалиться в штопор.
Симуляция подъёмной силы Ньютона методом частиц на CUDA