Skip to content

Commit

Permalink
Checkpoint.
Browse files Browse the repository at this point in the history
  • Loading branch information
edisonlee0212 committed Oct 18, 2024
1 parent 1da5452 commit 2d6bacd
Show file tree
Hide file tree
Showing 6 changed files with 53 additions and 29 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -84,7 +84,7 @@ mat3 cross_product_matrix(in vec3 v) {
return mat3(0, -v[2], v[1], v[2], 0, -v[0], -v[1], v[0], 0);
}

void solve_ldlt(in mat6 A, out vec6 x, out vec6 y, out vec6 z, in vec6 b) {
void solve_ldlt(in mat6 A, out vec6 x, in vec6 b) {
mat6 L;
L.v[0].v[0] = 1.0;
L.v[0].v[1] = 0.0;
Expand Down Expand Up @@ -136,6 +136,7 @@ void solve_ldlt(in mat6 A, out vec6 x, out vec6 y, out vec6 z, in vec6 b) {
D.v[4] = 0.0;
D.v[5] = 0.0;

vec6 y, z;

// Decomposition
for (int i = 0; i < 6; i++) {
Expand Down Expand Up @@ -292,9 +293,9 @@ void solve_stretch_bending_twisting_constraints(in int constraint_handle, in flo
jmjt.v[4].v[4] += bending_and_torsion_compliance[1];
jmjt.v[5].v[5] += bending_and_torsion_compliance[2];

vec6 delta_lambda, y, z;
vec6 delta_lambda;

solve_ldlt(jmjt, delta_lambda, y, z, rhs);
solve_ldlt(jmjt, delta_lambda, rhs);

lambda_sum.v[0] += delta_lambda.v[0];
lambda_sum.v[1] += delta_lambda.v[1];
Expand All @@ -315,21 +316,21 @@ void solve_stretch_bending_twisting_constraints(in int constraint_handle, in flo
x0_correction += inv_mass0 * delta_lambda_stretch;
vec3 v =
inverse_inertia0 * r0_cross_t * (-1.0 * delta_lambda_stretch) + m_inv_jt0 * delta_lambda_bending_and_torsion;
q0_correction[0] = g0[0][0] * v[0] + g0[0][1] * v[1] + g0[0][2] * v[2];
q0_correction[1] = g0[1][0] * v[0] + g0[1][1] * v[1] + g0[1][2] * v[2];
q0_correction[2] = g0[2][0] * v[0] + g0[2][1] * v[1] + g0[2][2] * v[2];
q0_correction[3] = g0[3][0] * v[0] + g0[3][1] * v[1] + g0[3][2] * v[2];
//q0_correction[0] = g0[0][0] * v[0] + g0[0][1] * v[1] + g0[0][2] * v[2];
//q0_correction[1] = g0[1][0] * v[0] + g0[1][1] * v[1] + g0[1][2] * v[2];
//q0_correction[2] = g0[2][0] * v[0] + g0[2][1] * v[1] + g0[2][2] * v[2];
//q0_correction[3] = g0[3][0] * v[0] + g0[3][1] * v[1] + g0[3][2] * v[2];
}

if (inv_mass1 != 0.0) {
x1_correction -= inv_mass1 * delta_lambda_stretch;
vec3 v = inverse_inertia1 * r1_cross_t * delta_lambda_stretch + m_inv_jt1 * delta_lambda_bending_and_torsion;
q1_correction[0] = g1[0][0] * v[0] + g1[0][1] * v[1] + g1[0][2] * v[2];
q1_correction[1] = g1[1][0] * v[0] + g1[1][1] * v[1] + g1[1][2] * v[2];
q1_correction[2] = g1[2][0] * v[0] + g1[2][1] * v[1] + g1[2][2] * v[2];
q1_correction[3] = g1[3][0] * v[0] + g1[3][1] * v[1] + g1[3][2] * v[2];
//q1_correction[0] = g1[0][0] * v[0] + g1[0][1] * v[1] + g1[0][2] * v[2];
//q1_correction[1] = g1[1][0] * v[0] + g1[1][1] * v[1] + g1[1][2] * v[2];
//q1_correction[2] = g1[2][0] * v[0] + g1[2][1] * v[1] + g1[2][2] * v[2];
//q1_correction[3] = g1[3][0] * v[0] + g1[3][1] * v[1] + g1[3][2] * v[2];
}

/*

RodConstraint rod_constraint = rod_constraints[constraint_handle];
rod_constraint.omega.xyz = omega;
Expand Down Expand Up @@ -376,14 +377,23 @@ void solve_stretch_bending_twisting_constraints(in int constraint_handle, in flo
rod_constraint.z[4] = z.v[4];
rod_constraint.z[5] = z.v[5];

rod_constraint.D[0] = D.v[0];
rod_constraint.D[1] = D.v[1];
rod_constraint.D[2] = D.v[2];
rod_constraint.D[3] = D.v[3];
rod_constraint.D[4] = D.v[4];
rod_constraint.D[5] = D.v[5];


rod_constraint.rhs[0] = rhs.v[0];
rod_constraint.rhs[1] = rhs.v[1];
rod_constraint.rhs[2] = rhs.v[2];
rod_constraint.rhs[3] = rhs.v[3];
rod_constraint.rhs[4] = rhs.v[4];
rod_constraint.rhs[5] = rhs.v[5];

rod_constraint.A = jmjt;
rod_constraints[constraint_handle] = rod_constraint;
*/
}


Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,14 @@ struct RodProperties {
float padding2;
};

struct vec6 {
float v[6];
};

struct mat6 {
vec6 v[6];
};

layout(set = DYNAMIC_STRANDS_PHYSICS_SET, binding = 0) readonly uniform ROD_PROPERTIES_BLOCK {
RodProperties rod_properties;
};
Expand Down Expand Up @@ -52,6 +60,7 @@ struct RodConstraint {
float padding0;
float padding1;

/*
vec4 omega;
mat4 j_omega0;
Expand All @@ -78,10 +87,13 @@ struct RodConstraint {
vec4 q0_correction;
vec4 q1_correction;
mat6 A;
float x[8];
float y[8];
float z[8];
float D[8];
float rhs[8];
*/
};

layout(std430, set = DYNAMIC_STRANDS_PHYSICS_SET, binding = 2) buffer ROD_CONSTRAINT_BLOCK {
Expand All @@ -96,13 +108,7 @@ layout(std430, set = DYNAMIC_STRANDS_PHYSICS_SET, binding = 3) buffer PER_SEGMEN
PerSegmentData per_segment_data_list[];
};

struct vec6 {
float v[6];
};

struct mat6 {
vec6 v[6];
};


mat3 compute_inverse_inertia_tensor_w(in StrandSegment strand_segment) {
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -192,13 +192,20 @@ void DynamicStrands::InitializeStrandsGroup(const InitializeParameters& initiali
segment.mass = segment.radius * segment.radius * glm::pi<float>() * initialize_parameters.wood_density;
segment.inv_mass = 1.f / segment.mass;
segment.inertia_tensor = ComputeInertiaTensorRod(segment.mass, segment.radius, segment.length);
segment.inertia_tensor = glm::vec3(1);

segment.inv_inertia_tensor = 1.f / segment.inertia_tensor;
}

segment.length = strand_group.GetStrandSegmentLength(static_cast<int>(i));
});
Upload();

const auto frame_count = Platform::GetMaxFramesInFlight();
for (uint32_t current_frame_index = 0; current_frame_index < frame_count; current_frame_index++) {
strands_descriptor_sets[current_frame_index]->UpdateBufferDescriptorBinding(0, device_strand_segments_buffer, 0);
strands_descriptor_sets[current_frame_index]->UpdateBufferDescriptorBinding(1, device_strands_buffer, 0);
}
InitializeOperators(initialize_parameters);
InitializeConstraints(initialize_parameters);
}
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -160,7 +160,7 @@ class DsStiffRod final : public IDynamicStrandsConstraint {
//---------
//| Debug |
//---------

/*
glm::vec3 omega;
float padding_omega;
Expand Down Expand Up @@ -195,12 +195,13 @@ class DsStiffRod final : public IDynamicStrandsConstraint {
glm::quat q0_correction;
glm::quat q1_correction;

mat6 A;
float x[8];
float y[8];
float z[8];
float D[8];
float rhs[8];

*/
};

struct PerStrandData {
Expand Down
5 changes: 0 additions & 5 deletions EvoEngine_Plugins/EcoSysLab/src/DynamicStrands.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -37,11 +37,6 @@ DynamicStrands::DynamicStrands() {
void DynamicStrands::Step(const StepParameters& target_step_parameters) {
if (strand_segments.empty())
return;
const auto current_frame_index = Platform::GetCurrentFrameIndex();

strands_descriptor_sets[current_frame_index]->UpdateBufferDescriptorBinding(0, device_strand_segments_buffer, 0);
strands_descriptor_sets[current_frame_index]->UpdateBufferDescriptorBinding(1, device_strands_buffer, 0);

if (target_step_parameters.physics) {
Physics(target_step_parameters.physics_parameters, operators, pre_step, constraints, post_step);
}
Expand Down
9 changes: 7 additions & 2 deletions EvoEngine_Plugins/EcoSysLab/src/DynamicStrandsPhysics.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -477,6 +477,8 @@ void DsStiffRod::InitializeData(const DynamicStrands::InitializeParameters& init
const float segment_length = (l0 + l1) * .5f;
rod_constraint.rest_darboux_vector = ComputeDarbouxVector(q0, q1, segment_length);
rod_constraint.average_segment_length = segment_length;


const auto rotation0_transpose = glm::transpose(glm::mat3_cast(q0));
const auto rotation1_transpose = glm::transpose(glm::mat3_cast(q1));
const auto& p0 = segment0.x0;
Expand All @@ -500,7 +502,7 @@ void DsStiffRod::InitializeData(const DynamicStrands::InitializeParameters& init
const auto torsional_stiffness = shear_modulus * polar_moment_of_inertia / segment_length;

rod_constraint.stiffness_coefficient_k = glm::vec3(stretching_stiffness, bending_stiffness, torsional_stiffness);

rod_constraint.stiffness_coefficient_k = glm::vec3(0.9, 0.9, 0.9);
segment_handle = next_segment_handle;
next_segment_handle = target_dynamic_strands.strand_segments[segment_handle].next_handle;
}
Expand Down Expand Up @@ -1004,6 +1006,9 @@ void DsStiffRod::Project(const DynamicStrands::PhysicsParameters& physics_parame
q1_correction.w = g1[3][0] * v[0] + g1[3][1] * v[1] + g1[3][2] * v[2];
}


//Debug
/*
GpuRodConstraint rod_constraint = rod_constraints[constraint_handle];
rod_constraint.omega = omega;
rod_constraint.j_omega0 = j_omega0;
Expand Down Expand Up @@ -1042,7 +1047,7 @@ void DsStiffRod::Project(const DynamicStrands::PhysicsParameters& physics_parame
rod_constraint.rhs[4] = rhs.v[4];
rod_constraint.rhs[5] = rhs.v[5];
rod_constraints[constraint_handle] = rod_constraint;
rod_constraints[constraint_handle] = rod_constraint;*/
};

const auto project_constraint = [&](const int rod_constraint_handle) {
Expand Down

0 comments on commit 2d6bacd

Please sign in to comment.