Skip to content

Commit

Permalink
multi-mesh VLM now working !
Browse files Browse the repository at this point in the history 
- added basic mesh files for multi mesh testing
- added basic multi mesh tests
- reenabled cd calculation
  • Loading branch information
@samayala22
samayala22 committed Aug 8, 2024
1 parent c6eb145 commit 5e395d6
Show file tree
Hide file tree
Showing 12 changed files with 169 additions and 124 deletions.
1 change: 1 addition & 0 deletions mesh/infinite_rectangular_2x2.x
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@
0.0000000000000000 0.0000000000000000 0.0000000000000000 20000.000000000000
20000.000000000000 20000.000000000000 40000.000000000000 40000.000000000000
40000.000000000000

0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000
13 changes: 13 additions & 0 deletions mesh/infinite_rectangular_2x2_part0.x
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,13 @@
1
3 3 1
0.0000000000000000 0.50000000000000000 1.0000000000000000
0.0000000000000000 0.50000000000000000 1.0000000000000000
0.0000000000000000 0.50000000000000000 1.0000000000000000

0.0000000000000000 0.00000000000000000 0.0000000000000000
10000.000000000000 10000.000000000000 10000.000000000000
20000.000000000000 20000.000000000000 20000.000000000000

0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.0000000000000000 0.0000000000000000
13 changes: 13 additions & 0 deletions mesh/infinite_rectangular_2x2_part1.x
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,13 @@
1
3 3 1
0.0000000000000000 0.50000000000000000 1.0000000000000000
0.0000000000000000 0.50000000000000000 1.0000000000000000
0.0000000000000000 0.50000000000000000 1.0000000000000000

20000.000000000000 20000.000000000000 20000.000000000000
30000.000000000000 30000.000000000000 30000.000000000000
40000.000000000000 40000.000000000000 40000.000000000000

0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.0000000000000000 0.0000000000000000
19 changes: 19 additions & 0 deletions mesh/infinite_rectangular_2x4.x
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,19 @@
1
3 5 1
0.0000000000000000 0.50000000000000000 1.0000000000000000
0.0000000000000000 0.50000000000000000 1.0000000000000000
0.0000000000000000 0.50000000000000000 1.0000000000000000
0.0000000000000000 0.50000000000000000 1.0000000000000000
0.0000000000000000 0.50000000000000000 1.0000000000000000

0.0000000000000000 0.00000000000000000 0.0000000000000000
10000.000000000000 10000.000000000000 10000.000000000000
20000.000000000000 20000.000000000000 20000.000000000000
30000.000000000000 30000.000000000000 30000.000000000000
40000.000000000000 40000.000000000000 40000.000000000000

0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.0000000000000000 0.0000000000000000
72 changes: 1 addition & 71 deletions tests/test_vlm_elliptic_coeffs.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35,83 +35,13 @@ float compute_analytical_gamma(float y, float a, float b, float alpha) {
return gamma0 * std::sqrt(1.0f - ratio*ratio);
}

// bool elliptic_convergence() {
// tiny::Config cfg("../../../../config/elliptic.vlm");

// vlm::VLM vlm(cfg);

// std::vector<int> dimensions = {
// 16, 32, 45, 64, 90, 128
// };

// const float a = 1.0f; // wing chord root
// const float b = 5.0f; // half wing span

// std::vector<double> norm_l1;
// std::vector<double> norm_l2;
// std::vector<double> norm_linf;

// const float alpha = 0.1f; // degrees

// for (const auto& dim : dimensions) {
// std::string filename = std::format("../../../../mesh/elliptic_{}x{}.x", dim, dim);
// vlm.mesh.io_read(filename);
// vlm.init();
// vlm::Solver solver(vlm.mesh, vlm.data, cfg);
// solver.run(alpha);

// double l1 = 0.0f;
// double l2 = 0.0f;
// double linf = 0.0f;
// int begin = (vlm.mesh.nc - 1) * vlm.mesh.ns;
// int end = vlm.mesh.nc * vlm.mesh.ns;
// // loop over last row of panels
// for (int i = begin; i < end; i++) {
// const float y = vlm.mesh.colloc.y[i];
// const float gamma = vlm.data.gamma[i];
// const float gamma_analytical = analytical_gamma(y, a, b, alpha);
// const double error = std::abs((gamma - gamma_analytical) / (gamma_analytical + 1e-7f));
// std::printf("y: %f, gamma: %f, gamma_analytical: %f, error: %f \n", y, gamma, gamma_analytical, error);

// l1 += error;
// l2 += error * error;
// linf = std::max(linf, error);
// }
// l1 /= (end - begin);
// l2 = std::sqrt(l2) / (end - begin);
// std::printf("L1: %f, L2: %f, Linf: %f\n", l1, l2, linf);

// norm_l1.push_back(l1);
// norm_l2.push_back(l2);
// norm_linf.push_back(linf);
// }

// double order_l1 = 0.0f;
// double order_l2 = 0.0f;
// double order_linf = 0.0f;

// auto order = [=](double norm0, double norm1, float dim0, float dim1) {
// return std::log(norm0 / norm1) / std::log((b/dim0)/(b/dim1));
// };

// for (int i = 0; i < dimensions.size() - 1; i++) {
// order_l1 += order(norm_l1[i], norm_l1[i+1], dimensions[i], dimensions[i+1]);
// order_l2 += order(norm_l2[i], norm_l2[i+1], dimensions[i], dimensions[i+1]);
// order_linf += order(norm_linf[i], norm_linf[i+1], dimensions[i], dimensions[i+1]);
// }
// order_l1 /= (dimensions.size() - 1);
// order_l2 /= (dimensions.size() - 1);
// order_linf /= (dimensions.size() - 1);
// std::printf("Order L1: %f, Order L2: %f, Order Linf: %f\n", order_l1, order_l2, order_linf);
// return 0;
// }

int main(int /*argc*/, char ** /*argv*/) {
const float a = 1.0f; // wing chord root
const float b = 5.0f; // half wing span

const std::vector<std::string> meshes = {"../../../../mesh/elliptic_64x64.x"};
const std::vector<std::string> backends = get_available_backends();

std::vector<f32> test_alphas = {0, 1, 2, 3, 4, 5, 10, 15};
std::transform(test_alphas.begin(), test_alphas.end(), test_alphas.begin(), to_radians);

Expand Down
49 changes: 49 additions & 0 deletions tests/test_vlm_square_coeffs.cpp
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,49 @@
#include <vector>
#include <string>
#include <fstream>
#include <cmath>
#include <cstdio>

#include "tinycombination.hpp"

#include "vlm.hpp"
#include "vlm_utils.hpp"

using namespace vlm;

float compute_analytical_cl(float alpha) {
return 2.0f * vlm::PI_f * alpha;
}

int main(int /*argc*/, char ** /*argv*/) {
const std::vector<std::vector<std::string>> meshes = {
{
"../../../../mesh/infinite_rectangular_2x2_part0.x",
"../../../../mesh/infinite_rectangular_2x2_part1.x"
}
};

const std::vector<std::string> backends = get_available_backends();
std::vector<f32> test_alphas = {0, 1, 2, 3, 4, 5, 10, 15};
//std::vector<f32> test_alphas = {1.f};
std::transform(test_alphas.begin(), test_alphas.end(), test_alphas.begin(), to_radians);

auto solvers = tiny::make_combination(meshes, backends);
for (const auto& [meshes_names, backend_name] : solvers) {
std::printf(">>> BACKEND: %s\n", backend_name.get().c_str());

VLM simulation{backend_name, meshes_names};

for (u64 i = 0; i < test_alphas.size(); i++) {
const FlowData flow{test_alphas[i], 0.0f, 1.0f, 1.0f};
const auto coeffs = simulation.run(flow);
const f32 analytical_cl = compute_analytical_cl(flow.alpha);
const f32 cl_aerr = std::abs(coeffs.cl - analytical_cl);
const f32 cl_rerr = analytical_cl < EPS_f ? 0.f : cl_aerr / analytical_cl;
std::printf(">>> Alpha: %.1f | CL = %.7f CD = %.7f CMx = %.6f CMy = %.6f CMz = %.6f\n", to_degrees(flow.alpha), coeffs.cl, coeffs.cd, coeffs.cm.x, coeffs.cm.y, coeffs.cm.z);
std::printf(">>> Analytical CL: %.7f | Abs Error: %.3E | Relative Error: %.5f%% \n", analytical_cl, cl_aerr, cl_rerr*100.f);
if (cl_rerr > 0.01f) return 1;
}
}
return 0;
}
9 changes: 5 additions & 4 deletions vlm/backends/cpu/include/vlm_backend_cpu.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -9,11 +9,10 @@ class BackendCPU final : public Backend {
public:
BackendCPU();
~BackendCPU();
void lhs_assemble(View<f32, Matrix<MatrixLayout::ColMajor>>& lhs, const View<f32, MultiSurface>& colloc, const View<f32, MultiSurface>& normals, const View<f32, MultiSurface>& verts_wing, const View<f32, MultiSurface>& verts_wake, u32 iteration) override;
void lhs_assemble(View<f32, Matrix<MatrixLayout::ColMajor>>& lhs, const View<f32, MultiSurface>& colloc, const View<f32, MultiSurface>& normals, const View<f32, MultiSurface>& verts_wing, const View<f32, MultiSurface>& verts_wake, std::vector<u32>& condition, u32 iteration) override;
void rhs_assemble_velocities(View<f32, MultiSurface>& rhs, const View<f32, MultiSurface>& normals, const View<f32, MultiSurface>& velocities) override;
void rhs_assemble_wake_influence(View<f32, MultiSurface>& rhs, const View<f32, MultiSurface>& gamma, const View<f32, MultiSurface>& colloc, const View<f32, MultiSurface>& normals, const View<f32, MultiSurface>& verts_wake, u32 iteration) override;
void displace_wake_rollup(View<f32, MultiSurface>& wake_rollup, const View<f32, MultiSurface>& verts_wake, const View<f32, MultiSurface>& verts_wing, const View<f32, MultiSurface>& gamma_wing, const View<f32, MultiSurface>& gamma_wake, f32 dt, u32 iteration) override;
// void displace_wing(const linalg::alias::float4x4& transform) override;
void displace_wing(const std::vector<linalg::alias::float4x4>& transforms, View<f32, MultiSurface>& verts_wing, View<f32, MultiSurface>& verts_wing_init) override;
void wake_shed(const View<f32, MultiSurface>& verts_wing, View<f32, MultiSurface>& verts_wake, u32 iteration) override;
void gamma_shed(View<f32, MultiSurface>& gamma_wing, View<f32, MultiSurface>& gamma_wing_prev, View<f32, MultiSurface>& gamma_wake, u32 iteration) override;
Expand All @@ -24,11 +23,13 @@ class BackendCPU final : public Backend {

// Per mesh kernels
f32 coeff_steady_cl_single(const View<f32, SingleSurface>& verts_wing, const View<f32, SingleSurface>& gamma_delta, const FlowData& flow, f32 area) override;
f32 coeff_steady_cl_multi(const View<f32, MultiSurface>& verts_wing, const View<f32, MultiSurface>& areas, const View<f32, MultiSurface>& gamma_delta, const FlowData& flow) override;
f32 coeff_steady_cl_multi(const View<f32, MultiSurface>& verts_wing, const View<f32, MultiSurface>& gamma_delta, const FlowData& flow, const View<f32, MultiSurface>& areas) override;
f32 coeff_steady_cd_single(const View<f32, SingleSurface>& verts_wake, const View<f32, SingleSurface>& gamma_wake, const FlowData& flow, f32 area) override;
f32 coeff_steady_cd_multi(const View<f32, MultiSurface>& verts_wake, const View<f32, MultiSurface>& gamma_wake, const FlowData& flow, const View<f32, MultiSurface>& areas) override;

// f32 coeff_unsteady_cl_single(const View<f32, SingleSurface>& verts_wing, const View<f32, SingleSurface>& gamma_delta, const View<f32, SingleSurface>& gamma, const View<f32, SingleSurface>& gamma_prev, const View<f32, SingleSurface>& local_velocities, const View<f32, SingleSurface>& areas, const View<f32, SingleSurface>& normals, const linalg::alias::float3& freestream, f32 dt, f32 area) override;
// f32 coeff_unsteady_cl_multi(const View<f32, MultiSurface>& verts_wing, const View<f32, MultiSurface>& areas, const View<f32, MultiSurface>& gamma_delta, const View<f32, MultiSurface>& gamma, const View<f32, MultiSurface>& gamma_prev, const View<f32, MultiSurface>& local_velocities, const View<f32, MultiSurface>& normals, const linalg::alias::float3& freestream, f32 dt) override;
// linalg::alias::float3 coeff_steady_cm(const FlowData& flow, const f32 area, const f32 chord, const u64 j, const u64 n) override;
// f32 coeff_steady_cd(const FlowData& flow, const f32 area, const u64 j, const u64 n) override;

void mesh_metrics(const f32 alpha_rad, const View<f32, MultiSurface>& verts_wing, View<f32, MultiSurface>& colloc, View<f32, MultiSurface>& normals, View<f32, MultiSurface>& areas) override;
f32 mesh_mac(const View<f32, SingleSurface>& verts_wing, const View<f32, SingleSurface>& areas) override;
Expand Down
Loading

0 comments on commit 5e395d6

Please sign in to comment.