commit 0771562bf2d4cf3a36c21260b3bc23bc02a7d95d
parent 591a2919fcddf86dc8e2073f492658d2a54a19d2
Author: Julio C. Ramirez-Ceballos <julio@jcrcx.xyz>
Date: Thu, 28 Nov 2024 08:12:20 +0100
Updates
Diffstat:
6 files changed, 100 insertions(+), 107 deletions(-)
diff --git a/.editorconfig b/.editorconfig
diff --git a/clt.h b/clt.h
@@ -35,16 +35,15 @@ typedef struct {
/*void calc_stressCLT(MaterialProperties *mat_list, Laminate *lam, */
/* double F[6], int *fail_list, double dT, double dM);*/
-void assemble_Z(double *Z, Laminate *lam);
+void assemble_Z(double *Z, Laminate *lam, bool debug);
void calc_thermal_forces(double *Nt, MaterialProperties *mat_list,
Laminate *lam, double *Z, int *fail_list, double dT);
void calc_moisture_forces(double *Nm, MaterialProperties *mat_list,
Laminate *lam, double *Z, int *fail_list, double dM);
-void assemble_Q(double Q[3][3], MaterialProperties *mat_prop,
- int fail_type);
+void assemble_Q(double Q[3][3], MaterialProperties *mat_prop);
void assemble_T(double T[3][3], double angle);
void assemble_ABD(double ABD[6][6], MaterialProperties *mat_list,
- Laminate *lam, double *Z, int *fail_list);
+ Laminate *lam, double *Z);
void fs_hashin_2D(MaterialProperties *mat_prop, double sig1, double sig2,
double tau, FSResult *result);
void hashin_2D(MaterialProperties *mat_list, Laminate *lam,
@@ -64,7 +63,7 @@ void maxstress_2D(MaterialProperties *mat_list, Laminate *lam,
void print_material_properties(MaterialProperties *mat, int mat_id);
void
-assemble_Z(double *Z, Laminate *lam)
+assemble_Z(double *Z, Laminate *lam, bool debug)
{
int i;
double total_thk = 0.0;
@@ -77,6 +76,13 @@ assemble_Z(double *Z, Laminate *lam)
for (i = 1; i <= lam->num_layers; ++i) {
Z[i] = Z[i-1] + lam->thk[i-1];
}
+
+ if (debug) {
+ printf("Layer Interface Positions (Z):\n");
+ for (i = 0; i <= lam->num_layers; ++i) {
+ printf("Z[%d] = %.4e mm\n", i, Z[i]);
+ }
+ }
}
void
@@ -101,7 +107,7 @@ calc_thermal_forces(double *Nt, MaterialProperties *mat_list,
}
/* Assemble Q matrix */
- assemble_Q(Q, mat_prop, fail_list ? fail_list[i] : 0);
+ assemble_Q(Q, mat_prop);
/* Calculate Qbar */
transpose(T,Ti);
@@ -143,7 +149,7 @@ calc_moisture_forces(double *Nm, MaterialProperties *mat_list,
}
/* Assemble Q matrix */
- assemble_Q(Q, mat_prop, fail_list ? fail_list[i] : 0);
+ assemble_Q(Q, mat_prop);
/* Calculate Qbar */
transpose(T,Ti);
@@ -163,7 +169,7 @@ calc_moisture_forces(double *Nm, MaterialProperties *mat_list,
}
void
-assemble_Q(double Q[3][3], MaterialProperties *mat_prop, int fail_type)
+assemble_Q(double Q[3][3], MaterialProperties *mat_prop)
{
if (mat_prop->isotropic) {
@@ -171,12 +177,12 @@ assemble_Q(double Q[3][3], MaterialProperties *mat_prop, int fail_type)
E = mat_prop->E1;
nu = mat_prop->nu12;
- if (fail_type == 1 || fail_type == 2 || fail_type == 3) {
- double df = 0.001;
- E *= df;
- nu *= df;
- }
-
+ /*if (fail_type == 1 || fail_type == 2 || fail_type == 3) {*/
+ /* double df = 0.001; */
+ /* E *= df;*/
+ /* nu *= df;*/
+ /*}*/
+ /**/
E_div = E / (1 - nu * nu);
Q[0][0] = Q[1][1] = E_div;
@@ -190,16 +196,18 @@ assemble_Q(double Q[3][3], MaterialProperties *mat_prop, int fail_type)
G12 = mat_prop->G12;
double df = 0.001; /* Degradation factor */
- if (fail_type == 1 || fail_type == 2) { /* fiber or shear failure */
- E1 *= df;
- E2 *= df;
- nu12 *= df;
- G12 *= df;
- } else if (fail_type == 3) { /* matrix failure */
- E2 *= df;
- nu12 *= df;
- G12 *= df;
- }
+ /* fiber or shear failure */
+ /*if (fail_type == 1 || fail_type == 2) { */
+ /* E1 *= df; */
+ /* E2 *= df; */
+ /* nu12 *= df; */
+ /* G12 *= df;*/
+ /*} else if (fail_type == 3) { */
+ /* matrix failure */
+ /* E2 *= df; */
+ /* nu12 *= df; */
+ /* G12 *= df;*/
+ /*}*/
double n21 = nu12 * E2 / E1;
@@ -226,7 +234,7 @@ assemble_T(double T[3][3], double angle)
void
assemble_ABD(double ABD[6][6], MaterialProperties *mat_list,
- Laminate *lam, double *Z, int *fail_list)
+ Laminate *lam, double *Z)
{
double A[3][3] = {0}, B[3][3] = {0}, D[3][3] = {0};
double Q[3][3], T[3][3], Ti[3][3];
@@ -235,7 +243,7 @@ assemble_ABD(double ABD[6][6], MaterialProperties *mat_list,
for (int i = 0; i < lam->num_layers; ++i) {
MaterialProperties *mat_prop = &mat_list[lam->mat_id[i]];
- assemble_Q(Q, mat_prop, fail_list ? fail_list[i] : 0);
+ assemble_Q(Q, mat_prop);
assemble_T(T, lam->ang[i]);
/* Transpose T to get Ti */
@@ -599,35 +607,31 @@ double nu_xy(double ABD[6][6], double h) {
double matrix1[5][5], matrix2[5][5];
double det1, det2, nxy;
- matrix1[0][0] = ABD[1][0];
- matrix1[0][1] = ABD[1][2];
- matrix1[0][2] = ABD[1][3];
- matrix1[0][3] = ABD[1][4];
- matrix1[0][4] = ABD[1][5];
-
- matrix1[1][0] = ABD[2][0];
- matrix1[1][1] = ABD[2][2];
- matrix1[1][2] = ABD[2][3];
- matrix1[1][3] = ABD[2][4];
- matrix1[1][4] = ABD[2][5];
-
- matrix1[2][0] = ABD[3][0];
- matrix1[2][1] = ABD[3][2];
- matrix1[2][2] = ABD[3][3];
- matrix1[2][3] = ABD[3][4];
- matrix1[2][4] = ABD[3][5];
-
- matrix1[3][0] = ABD[4][0];
- matrix1[3][1] = ABD[4][2];
- matrix1[3][2] = ABD[4][3];
- matrix1[3][3] = ABD[4][4];
- matrix1[3][4] = ABD[4][5];
-
- matrix1[4][0] = ABD[5][0];
- matrix1[4][1] = ABD[5][2];
- matrix1[4][2] = ABD[5][3];
- matrix1[4][3] = ABD[5][4];
- matrix1[4][4] = ABD[5][5];
+ matrix1[0][0] = ABD[1][0];
+ matrix1[0][1] = ABD[1][2];
+ matrix1[0][2] = ABD[1][3];
+ matrix1[0][3] = ABD[1][4];
+ matrix1[0][4] = ABD[1][5];
+ matrix1[1][0] = ABD[2][0];
+ matrix1[1][1] = ABD[2][2];
+ matrix1[1][2] = ABD[2][3];
+ matrix1[1][3] = ABD[2][4];
+ matrix1[1][4] = ABD[2][5];
+ matrix1[2][0] = ABD[3][0];
+ matrix1[2][1] = ABD[3][2];
+ matrix1[2][2] = ABD[3][3];
+ matrix1[2][3] = ABD[3][4];
+ matrix1[2][4] = ABD[3][5];
+ matrix1[3][0] = ABD[4][0];
+ matrix1[3][1] = ABD[4][2];
+ matrix1[3][2] = ABD[4][3];
+ matrix1[3][3] = ABD[4][4];
+ matrix1[3][4] = ABD[4][5];
+ matrix1[4][0] = ABD[5][0];
+ matrix1[4][1] = ABD[5][2];
+ matrix1[4][2] = ABD[5][3];
+ matrix1[4][3] = ABD[5][4];
+ matrix1[4][4] = ABD[5][5];
/* Perform LU decomposition for matrix1 */
if (LUDecompose(&matrix1[0][0], 5, 5) < 0) {
@@ -637,6 +641,7 @@ double nu_xy(double ABD[6][6], double h) {
/* Calculate determinant for matrix1 */
det1 = LUDeterminant(&matrix1[0][0], 5, 5);
+ printf("%.4e\n",det1);
matrix2[0][0] = ABD[1][1];
matrix2[0][1] = ABD[1][2];
@@ -666,21 +671,22 @@ double nu_xy(double ABD[6][6], double h) {
// Perform LU decomposition for matrix2
if (LUDecompose(&matrix2[0][0], 5, 5) < 0) {
- fprintf(stderr, "Error during LU decomposition for matrix2.\n");
+ fprintf(stderr, "Error during LU decomposition for submatrix.\n");
return -1;
}
// Calculate determinant for matrix2
det2 = LUDeterminant(&matrix2[0][0], 5, 5);
+ printf("%.4e\n",det2);
// Check for division by zero
if (det2 == 0) {
- fprintf(stderr, "Error: Division by zero in E_x calculation, determinant of the submatrix is zero.\n");
+ fprintf(stderr, "Error: Division by zero in nu_xy calculation, determinant of the submatrix is zero.\n");
return -1;
}
// Compute E_x
- nxy = -det1 / det2;
+ nxy = det1 / det2;
return nxy;
}
double nu_yx(double ABD[6][6], double h) {
@@ -692,21 +698,25 @@ double nu_yx(double ABD[6][6], double h) {
matrix1[0][2] = ABD[0][3];
matrix1[0][3] = ABD[0][4];
matrix1[0][4] = ABD[0][5];
+
matrix1[1][0] = ABD[2][0];
matrix1[1][1] = ABD[2][2];
matrix1[1][2] = ABD[2][3];
matrix1[1][3] = ABD[2][4];
matrix1[1][4] = ABD[2][5];
+
matrix1[2][0] = ABD[3][0];
matrix1[2][1] = ABD[3][2];
matrix1[2][2] = ABD[3][3];
matrix1[2][3] = ABD[3][4];
matrix1[2][4] = ABD[3][5];
+
matrix1[3][0] = ABD[4][0];
matrix1[3][1] = ABD[4][2];
matrix1[3][2] = ABD[4][3];
matrix1[3][3] = ABD[4][4];
matrix1[3][4] = ABD[4][5];
+
matrix1[4][0] = ABD[5][0];
matrix1[4][1] = ABD[5][2];
matrix1[4][2] = ABD[5][3];
@@ -759,12 +769,12 @@ double nu_yx(double ABD[6][6], double h) {
// Check for division by zero
if (det2 == 0) {
- fprintf(stderr, "Error: Division by zero in E_x calculation, determinant of the submatrix is zero.\n");
+ fprintf(stderr, "Error: Division by zero in E_y calculation, determinant of the submatrix is zero.\n");
return -1;
}
// Compute E_x
- nyx = -det1 / det2;
+ nyx = det1 / det2;
return nyx;
}
@@ -775,14 +785,20 @@ double nu_yx(double ABD[6][6], double h) {
/* StressStrainResult *result*/
/* ) */
/*{*/
-/* int num = lam->num_layers;*/
-/* double *Z = malloc((num + 1) * sizeof(double));*/
-/* assemble_Z(Z, &lam);*/
-/**/
-/* double h; for(i = 0; i < num; i++) h += lam.thk[i];*/
-/**/
-/* double ABD[N][N];*/
-/* assemble_ABD(ABD, materials, &lam, Z, NULL);*/
/**/
-/* printf("%d\n", NELEMS(Z));*/
/*}*/
+
+
+void print_equivalent_properties(double ABD[6][6], double h) {
+ double Ex = E_x(ABD, h);
+ double Ey = E_y(ABD, h);
+ double Gxy = G_xy(ABD, h);
+ double nxy = nu_xy(ABD, h);
+ double nyx = nu_yx(ABD, h);
+
+ printf("Equivalent elastic modulus Ex:\t%.4e\n", Ex);
+ printf("Equivalent elastic modulus Ey:\t%.4e\n", Ey);
+ printf("Equivalent shear modulus Gxy:\t%.4e\n", Gxy);
+ printf("Equivalent poisson ratio nu_xy:\t%.2f\n", nxy);
+ printf("Equivalent poisson ratio nu_yx:\t%.2f\n", nyx);
+}
diff --git a/laminate.dat b/laminate.dat
@@ -1,7 +1,3 @@
-1 0.127 0
-1 0.127 90
0 0.127 0
-2 0.127 0
0 0.127 0
-1 0.127 90
-1 0.127 0
+0 0.127 0
diff --git a/main b/main
Binary files differ.
diff --git a/main.c b/main.c
@@ -2,7 +2,7 @@
* main.c
*
* Author: Julio C. Ramirez-Ceballos
- * Licence: MIT
+ * Licence: ISC
* Description:
* This code is used to calculate the equivalent properties of a laminate
* with some configuration.
@@ -26,7 +26,7 @@ int main()
double F[] = {
100000,
0e2,
- 0e4,
+ 100e4,
10e0,
0e0,
0e0,
@@ -37,20 +37,20 @@ int main()
read_material_properties("materials.dat", &num_materials, true);
Laminate lam;
- read_laminate_data("laminate.dat", &lam, true);
+ read_laminate_data("laminate.dat", &lam, false);
double *Z = malloc((lam.num_layers + 1) * sizeof(double));
- assemble_Z(Z, &lam);
- printf("Layer Interface Positions (Z):\n");
- for (i = 0; i <= lam.num_layers; ++i) {
- printf("Z[%d] = %.4e mm\n", i, Z[i]);
- }
+ assemble_Z(Z, &lam, false);
printf("\n");
double ABD[N][N];
- assemble_ABD(ABD, materials, &lam, Z, NULL);
+ assemble_ABD(ABD, materials, &lam, Z);
+ double h; for(i = 0; i < lam.num_layers; i++) h += lam.thk[i];
+ printf("h = %0.4e mm\n", h);
+
+ print_equivalent_properties(ABD, h);
if (LUDecompose(&ABD[0][0], 6, 6) < 0) {
fprintf(stderr, "Error during LU decomposition for matrix1.\n");
@@ -108,7 +108,7 @@ int main()
/*CORRIGER*/
/*for (i = 0; i < lam.num_layers; i++) {*/
/* MaterialProperties *mat_prop = &materials[lam.mat_id[i]];*/
- /* assemble_Q(Q, mat_prop, fail_list ? fail_list[i] : 0);*/
+ /* assemble_Q(Q, mat_prop);*/
/* assemble_T(T, lam.ang[i]);*/
/* double temp[3];*/
/**/
@@ -126,27 +126,6 @@ int main()
/*}*/
/*printf("%.4e\n", MS_stress_inf[0][0]);*/
- double h;
- for(i = 0; i < lam.num_layers; i++) {
- h += lam.thk[i];
- }
- printf("h = %0.4e mm\n", h);
-
- double Ex = E_x(ABD, h);
- printf("Equivalent elastic modulus Ex:\t%.4e\n", Ex);
-
- double Ey = E_y(ABD, h);
- printf("Equivalent elastic modulus Ey:\t%.4e\n", Ey);
-
- double Gxy = G_xy(ABD, h);
- printf("Equivalent shear modulus Gxy:\t%.4e\n", Gxy);
-
- double nxy = nu_xy(ABD, h);
- printf("Equivalent poisson ratio nu_xy:\t%.2f\n", nxy);
-
- double nyx = nu_yx(ABD, h);
- printf("Equivalent poisson ratio nu_yx:\t%.2f\n", nyx);
-
free(lam.thk);
free(lam.ang);
free(lam.mat_id);
diff --git a/rw.c b/rw.c
@@ -68,7 +68,8 @@ void read_laminate_data(const char *filename, Laminate *lam, bool print_details)
double thickness, angle;
// Use sscanf to parse the line
- if (sscanf(line, "%d %lf %lf", &mat_id, &thickness, &angle) == 3) {
+ if (sscanf(line, "%d %lf %lf",
+ &mat_id, &thickness, &angle) == 3) {
lam->mat_id[i] = mat_id;
lam->thk[i] = thickness;
lam->ang[i] = angle;
@@ -111,12 +112,13 @@ read_material_properties(const char *filename,
int material_count = 0;
MaterialProperties *mat = NULL;
+ if (print_details){
printf("\n");
printf("+---------------------+\n");
printf("| Material Properties |\n");
printf("+---------------------+\n");
printf("\n");
-
+ }
file = fopen(filename, "r");
if (file == NULL) {
fprintf(stderr, "Error opening file %s\n", filename);
