yggdrasil75/stupidsimcpp
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

pushing this home. its working better now.

Browse Source
This commit is contained in:
Yggdrasil75
2025-11-28 13:13:05 -05:00
parent b429b1b98d
commit 439e0c0a2a
3 changed files with 55 additions and 55 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
tests/g2temp.cpp
View File

@@ -286,7 +286,7 @@ void mainLogic(const AnimationConfig& config, Shared& state, int gradnoise) {
}
//expandPixel(grid,config,seeds);
grid.diffuseTempsOptimized(100);
grid.diffuseTemps(100);
std::lock_guard<std::mutex> lock(state.mutex);
state.grid = grid;

64
util/grid/grid2.hpp
View File

@@ -977,25 +977,25 @@ public:
}
}
void diffuseTemps(int timestep) {
TIME_FUNCTION;
std::cout << "diffusing temps with a timestep of " << timestep << std::endl;
if (tempMap.empty() || timestep < 1) return;
#pragma omp parallel for
for (const auto& [id, pos] : Positions) {
auto tempIT = tempMap.find(id);
if (tempIT != tempMap.end()) {
Temp oldtemp = tempIT->second;
tempMap.erase(id);
float newtemp = Temp::calGrad(pos, getTemps(id));
float newtempMult = (newtemp-oldtemp.temp) * timestep;
oldtemp.temp = newtempMult;
tempMap.emplace(id, oldtemp);
}
}
}
// void diffuseTemps(int timestep) {
// TIME_FUNCTION;
// std::cout << "diffusing temps with a timestep of " << timestep << std::endl;
// if (tempMap.empty() || timestep < 1) return;
// #pragma omp parallel for
// for (const auto& [id, pos] : Positions) {
// auto tempIT = tempMap.find(id);
// if (tempIT != tempMap.end()) {
// Temp oldtemp = tempIT->second;
// tempMap.erase(id);
// float newtemp = Temp::calGrad(pos, getTemps(id));
// float newtempMult = (newtemp-oldtemp.temp) * timestep;
// oldtemp.temp = newtempMult;
// tempMap.emplace(id, oldtemp);
// }
// }
// }
void diffuseTempsOptimized(float deltaTime) {
void diffuseTemps(float deltaTime) {
TIME_FUNCTION;
if (tempMap.empty() || deltaTime <= 0) return;
@@ -1003,33 +1003,27 @@ public:
tempUpdates.reserve(tempMap.size());
// Process in spatial order for better cache performance
for (const auto& [id, tempObj] : tempMap) {
for (auto& [id, tempObj] : tempMap) {
Vec2 pos = Positions.at(id);
float oldtemp = tempObj.temp;
// Use smaller query radius for diffusion
auto nearbyIds = spatialGrid.queryRange(pos, neighborRadius * tempObj.conductivity);
float heatTransfer = 0.0f;
int validNeighbors = 0;
std::unordered_map<Vec2, Temp> neighborTemps;
for (size_t neighborId : nearbyIds) {
if (neighborId != id && tempMap.find(neighborId) != tempMap.end()) {
float tempDiff = tempMap.at(neighborId).temp - tempObj.temp;
float distance = pos.distance(Positions.at(neighborId));
if (distance > EPSILON) {
// Basic heat conduction formula
float conductivity = tempObj.conductivity * (tempObj.diffusivity + tempMap.at(neighborId).diffusivity);
heatTransfer += conductivity * tempDiff / distance;
validNeighbors++;
}
neighborTemps.emplace(Positions.at(neighborId), tempMap.at(neighborId));
}
}
if (validNeighbors > 0) {
float newTemp = tempObj.temp + heatTransfer * deltaTime / validNeighbors;
tempUpdates.emplace_back(id, newTemp);
}
tempObj.calGrad(pos, neighborTemps);
float newtemp = tempObj.temp;
float tempdiff = (oldtemp - newtemp) * (deltaTime / 1000);
// if (std::abs(newtemp) < EPSILON) {
// std::cout << "new temp is: " << newtemp << " old temp is: " << oldtemp << " diff*delta: " << tempdiff << std::endl;
// }
tempObj.temp = oldtemp - tempdiff;
}
// Batch update temperatures

42
util/simblocks/temp.hpp
View File

@@ -78,46 +78,52 @@ public:
return num / den;
}
static float calGrad(const Vec2& testPos, const std::unordered_map<Vec2, Temp>& others) {
void calGrad(const Vec2& testPos, const std::unordered_map<Vec2, Temp>& others) {
int meh;
//std::cout << meh++ << std::endl;
std::vector<std::pair<Vec2, float>> nearbyPoints;
for (const auto& [point, temp] : others) {
for (const auto& [point, ctemp] : others) {
if (point.distance(testPos) <= 25) {
nearbyPoints.emplace_back(point, temp.temp);
nearbyPoints.emplace_back(point, ctemp.temp);
}
}
float sumX, sumY, sumT, sumX2, sumY2, sumXY, sumXT, sumYT = 0;
//std::cout << meh++ << std::endl;
float sumX = 0,sumY = 0,sumT = 0,sumX2 = 0,sumY2 = 0,sumXY = 0,sumXT = 0, sumYT = 0;
int n = nearbyPoints.size();
for (const auto& [point, temp] : nearbyPoints) {
//std::cout << meh++ << std::endl;
for (const auto& [point, ctemp] : nearbyPoints) {
float x = point.x - testPos.x;
float y = point.y - testPos.y;
sumX += x;
sumY += y;
sumT += temp;
sumT += ctemp;
sumX2 += x * x;
sumY2 += y * y;
sumXY += x * y;
sumXT += x * temp;
sumYT += y * temp;
sumXT += x * ctemp;
sumYT += y * ctemp;
}
//std::cout << meh++ << std::endl;
float det = sumX2 * sumY2 - sumXY * sumXY;
Vec2 calpoint;
if (std::abs(det) < 1e-10) {
Vec2 calpoint = Vec2(0, 0); // Singular matrix, cannot solve
calpoint = Vec2(0, 0); // Singular matrix, cannot solve
}
float a = (sumXT * sumY2 - sumYT * sumXY) / det;
float b = (sumX2 * sumYT - sumXY * sumXT) / det;
Vec2 calpoint = Vec2(a, b); // ∇T = (∂T/∂x, ∂T/∂y)
else {
float a = (sumXT * sumY2 - sumYT * sumXY) / det;
float b = (sumX2 * sumYT - sumXY * sumXT) / det;
calpoint = Vec2(a, b); // ∇T = (∂T/∂x, ∂T/∂y)
}
//std::cout << meh++ << std::endl;
Vec2 closest = findClosestPoint(testPos, others);
Vec2 displacement = testPos - closest;
float refTemp = others.at(closest).temp;
float estimatedTemp = refTemp + (calpoint.x * displacement.x + calpoint.y * displacement.y);
return estimatedTemp;
float estimatedTemp = temp + (calpoint.x * displacement.x + calpoint.y * displacement.y);
temp = estimatedTemp;
}
};