/* nag_mesh_dim2_gen_delaunay (d06abc) Example Program.
*
* Copyright 2024 Numerical Algorithms Group.
*
* Mark 30.0, 2024.
*/
#include <math.h>
#include <nag.h>
#include <stdio.h>
int main(void) {
/* Scalars */
Integer exit_status = 0;
double c, d, dnvint, r, s, theta, theta_i, tmp;
Integer i, i1, ic, itrace, j, k, nedge, nelt, npropa, nrae;
Integer nv, nvb, nvint, nvmax, reftk, nearest, nv_near, nelt_near;
/* Arrays */
double t[2];
char pmesh[2];
double *coor = 0, *weight = 0;
Integer *conn = 0, *edge = 0;
/* Nag Types */
NagError fail;
INIT_FAIL(fail);
exit_status = 0;
printf("nag_mesh_dim2_gen_delaunay (d06abc) Example Program Results\n\n");
/* Skip heading in data file */
scanf("%*[^\n] ");
/* Reading of the geometry */
scanf("%" NAG_IFMT "%" NAG_IFMT "%" NAG_IFMT "%*[^\n] ", &nrae, &nvint,
&nvmax);
scanf(" ' %1s '%*[^\n]", pmesh);
nvb = nvint + 2 * nrae;
nedge = nvb;
if (nvb > nvmax) {
printf("Problem with the array dimensions\n");
printf(" nvb nvmax %6" NAG_IFMT "%6" NAG_IFMT "\n", nvb, nvmax);
printf(" Please increase the value of nvmax\n");
exit_status = -1;
goto END;
}
/* Allocate memory */
if (!(coor = NAG_ALLOC(2 * nvmax, double)) ||
!(weight = NAG_ALLOC(nvint, double)) ||
!(conn = NAG_ALLOC(3 * (2 * nvmax + 5), Integer)) ||
!(edge = NAG_ALLOC(3 * nedge, Integer))) {
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
/* Coordinates of the boundary mesh vertices and boundary edges */
/* Circular outer boundary, radius 3 and centre (1,0) */
theta = 2.0 * X01AAC / ((double)nvint);
r = 3.0;
t[0] = 1.0;
theta_i = 0.0;
for (i = 0; i < nvint; ++i) {
coor[2 * i] = r * cos(theta_i) + t[0];
coor[2 * i + 1] = r * sin(theta_i);
theta_i = theta_i + theta;
}
/* Read data for aerofoil RAE 2822 */
for (i = nvint; i < nvint + nrae; ++i) {
scanf("%" NAG_IFMT "%lf%lf%*[^\n] ", &i1, &coor[2 * i], &coor[2 * i + 1]);
}
/* Transform RAE 2822 for secondary foil */
theta = X01AAC / 12.0;
c = cos(theta);
s = sin(theta);
ic = 2 * (nvint + nrae);
i1 = 2 * nvint;
/* Copy and rotate coordinates by theta = pi/12 */
for (i = 0; i < 2 * nrae; ++i) {
coor[ic + i] = coor[i1 + i];
}
i1 = ic;
for (i = 0; i < nrae; ++i) {
tmp = coor[i1];
coor[i1] = s * coor[i1 + 1] + c * tmp;
coor[i1 + 1] = c * coor[i1 + 1] - s * tmp;
i1 = i1 + 2;
}
/* Reduce by 0.4 and translate to distance 0.25 from intercept at (0.75,0) */
d = 0.4;
t[0] = 0.75 + 0.25 * c;
t[1] = -0.25 * s;
for (i = 0; i < nrae; ++i) {
coor[ic + 2 * i] = d * coor[ic + 2 * i] + t[0];
coor[ic + 2 * i + 1] = d * coor[ic + 2 * i + 1] + t[1];
}
/* Boundary edges */
i1 = 0;
for (i = 0; i < nedge; ++i) {
edge[i1] = i + 1;
edge[i1 + 1] = i + 2;
edge[i1 + 2] = 0;
i1 = i1 + 3;
}
/* Tie up end of three boundary edges */
edge[3 * nvint - 2] = 1;
edge[3 * (nvint + nrae) - 2] = nvint + 1;
edge[3 * nedge - 2] = nvint + nrae + 1;
/* Initialize mesh control parameters */
itrace = 0;
/* Generation of interior vertices on the RAE airfoils wake */
dnvint = 2.5 / (double)(nvint + 1);
i1 = 2 * nvb;
for (i = 0; i < nvint; ++i) {
coor[i1] = (double)(i + 1) * dnvint + 1.38;
coor[i1 + 1] = -0.27 * coor[i1] + 0.2;
i1 = i1 + 2;
weight[i] = 0.01;
}
/* Loop on the propagation coef */
nearest = 250;
for (j = 0; j < 4; ++j) {
switch (j) {
case 0:
npropa = -5;
break;
case 1:
npropa = -1;
break;
case 2:
npropa = 1;
break;
default:
npropa = 5;
}
/* Call to the 2D Delaunay-Voronoi mesh generator */
/* nag_mesh_dim2_gen_delaunay (d06abc).
* Generates a two-dimensional mesh using a Delaunay-Voronoi
* process
*/
nag_mesh_dim2_gen_delaunay(nvb, nvint, nvmax, nedge, edge, &nv, &nelt, coor,
conn, weight, npropa, itrace, 0, &fail);
if (fail.code == NE_NOERROR) {
if (pmesh[0] == 'N') {
printf(" Mesh characteristics with npropa =%6" NAG_IFMT "\n", npropa);
nv_near = ((nv + nearest / 2) / nearest) * nearest;
nelt_near = ((nelt + nearest / 2) / nearest) * nearest;
printf(" nv =%10" NAG_IFMT " to the nearest %3" NAG_IFMT "\n",
nv_near, nearest);
printf(" nelt =%10" NAG_IFMT " to the nearest %3" NAG_IFMT "\n",
nelt_near, nearest);
} else if (pmesh[0] == 'Y') {
/* Output the mesh in a form suitable for printing */
printf(" %10" NAG_IFMT " %10" NAG_IFMT "\n", nv, nelt);
for (i = 0; i < nv; ++i) {
printf(" %15.6e %15.6e \n", coor[2 * i], coor[2 * i + 1]);
}
reftk = 0;
for (k = 0; k < nelt; ++k) {
printf(" %10" NAG_IFMT " %10" NAG_IFMT " %10" NAG_IFMT ""
" %10" NAG_IFMT "\n",
conn[3 * k], conn[3 * k + 1], conn[3 * k + 2], reftk);
}
} else {
printf("Problem with the printing option Y or N\n");
exit_status = -1;
goto END;
}
} else {
printf("Error from nag_mesh_dim2_gen_delaunay (d06abc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
}
END:
NAG_FREE(coor);
NAG_FREE(weight);
NAG_FREE(conn);
NAG_FREE(edge);
return exit_status;
}