Generating the Regular Axis from Irregular Column Grids through Genetic Algorithm

Featured Application: This study proposes an auto-generation method of the regular axis, applied in digital documentation of cultural heritage. The procedure is repeatable, the results show the displacement of columns with visual expression and numerical analysis. Abstract: Historic building information modeling (HBIM) provides an index frame for digital documentation of the cultural heritage, as a continuous process of reverse engineering. The index frame should be a regular model, with a clear comprehension of each component; consequently, associated knowledge could be live-updated with the investigation progress. Therefore, the method of data registration stresses its importance. The axis is fundamental for Chinese traditional architecture as the basis for positioning all components in wooden structures. However, displacement often happens. To correct the displacement while modeling, the hypothetical axis should be determined ﬁrst. This paper thus proposes a method of generating the regular axis from irregular column grids and aims to develop an automatic solution that is repeatable and transplantable. The ﬁnite element modeling (FEM) abstracts the actual problem to enable numerical calculation. Starting from a candidate solution, the genetic algorithm calculates a solution closest to the expectation from the possible solutions in several minutes. The standard deviation is used to measure the amount of displacement based on the hypothetical axis, which is expected to be minimum. This method is compatible with most kinds of input data, e.g., point cloud, excel data, 2D drawing, mesh model, etc., and applied to a World Heritage Site in Qufu (Shandong, China). The results show the displacement of columns with visual expression and numerical analysis and prove that the proposed method is repeatable and traceable and can easily be applied to different projects by changing input parameters.


Review of Related Works
During the digital documentation of cultural heritage, displacements are usually found in wooden structures [1][2][3][4], resulting in difficulty of positioning while modeling.
Different from modern architecture, in Chinese traditional architecture, the axis is designed for locating column heads, rather than column feet. Because eave columns are inclined and raised in an arithmetic sequence, from center to corner, the plane projection of centroids connection line of eave column heads is rectangular, while the eave column feet are not.
This manufacturing method is documented in "Ying Zao Fa Shi" [5][6][7][8][9] (which states building regulations of Song Dynasty, A.D. 1103), and it can be seen in many heritage sites in East Asia nowadays. However, this method is difficult to verify, as all the wooden structures have lost their original positions.
As a continuous process of reverse engineering, i.e., regularized rebuild or ideal reconstruction [10][11][12][13], an HBIM model should be an index frame [14,15], a prototype library [16][17][18], or an ideal model [19][20][21], of which components could be comprehended with architectural semantics, and knowledge could be associated through multidisciplinary collaboration. The traditional method of digital documentation in China accumulates a large number of 2D drawings, yet mistakes occur while modeling with reference to these data.
As a result, the data acquired from 2D drawings have to be registered, e.g., a hypothetical orthogonal grid system. Each axis is calculated through a certain algorithm so that the average distance between each "as-found" column centroid and "as-designed" column centroid is expected to be the smallest. Once an ideal grid system is established, the numerical analysis could be applied.
In recent years, the method of algorithm modeling (or generative modeling) with visual programming tools is developing increasingly [22][23][24][25][26][27][28][29][30][31]. This technique provides traceable solutions which document the whole modeling procedure, including input parameters, process algorithms, and resulting models. Based on this foundation, the research that considers grids as a whole of the interrelated system cloud be further developed.
To create a whole column grid system, finite element modeling (FEM) with the utilization of genetic algorithm is chosen as the method, as the multi-goals (result axes) are known quantities while the genomes (input parameters) are unknown quantities. The FEM transforms the problem of infinite degrees of freedom in a continuous domain into the problem of finite degrees of freedom in a discrete domain, enabling numerical calculation [32][33][34][35][36]. Furthermore, compiling algorithms for mathematical representation, the genetic algorithm is available to implement for optimizing FEM, with a Grasshopper plug-in, the Galapagos generic solvers [37]. The Rhinoceros thus is chosen as the software platform for the automatic solution through the genetic algorithm with FEM. The Galapagos as an evolutionary solver in Grasshopper has been widely used for parametric design and result optimization [38][39][40][41][42][43][44][45][46].

Purpose and Significance
This paper aims to propose a method of automatically generating whole axes in column grids at the same time with visible algorithm nodes and thousands of iterations, which would be less possible by manual operating.
Starting with a pre-given candidate solution and point cloud, the Galapagos module tries to move and rotate all of the axes with the expectation that the average distance between each ideal axis intersection point and as-found column head centroid is the smallest. When applying to different kinds of heritage sites, one only needs to change input parameters by adjusting the numeric slider or removing non-participate columns, etc., and then the procedure would finish the calculation within less than an hour.
The proposed method consists of several steps, as shown in the flowchart (Figure 1). Appl. Sci. 2022, 12, x FOR PEER REVIEW 3 of 17

Finite Element Modeling Using Genetic Algorithm
Take the Kuiwen Pavilion as an example. It is located in the Temple of Confucius, which has been part of the UNESCO World Heritage Site "Temple and Cemetery of Confucius and the Kong Family Mansion in Qufu" since 1994. The Kuiwen Pavilion was built in the Imperial style, which means that the structural members should ideally form column grids with rectangle spaces.
However, based on a plane section of column heads from the point cloud, an irregular grid system shows as-found column heads together with beams ( Figure 2).

Finite Element Modeling Using Genetic Algorithm
Take the Kuiwen Pavilion as an example. It is located in the Temple of Confucius, which has been part of the UNESCO World Heritage Site "Temple and Cemetery of Confucius and the Kong Family Mansion in Qufu" since 1994. The Kuiwen Pavilion was built in the Imperial style, which means that the structural members should ideally form column grids with rectangle spaces.
However, based on a plane section of column heads from the point cloud, an irregular grid system shows as-found column heads together with beams ( Figure 2).
Column heads serve as joints of beams, and together with beams, they form a spatial framework, therefore it can be assumed that displacements are evenly distributed from their original place, i.e., they meet Laplace-Gaussian distribution. To bring this problem to numerical calculation, the column head could be simplified to one circle, which could be deconstructed into two parameters: diameter and a centroid (point).
To define an axis, at least two centroids are required, yet this axis is certainly determined by axiom. However, if there are three points or more, an axis is uncertainly determined by least-squares fitting. Consequently, it can be assumed that grids are rectangles, displacements of column heads are evenly distributed along each rebuild axis, and the variance of distances between current centroids and rebuild centroids of column heads has a minimum value. Therefore, Gauss sum instead of variance is expected to be minimum; thus, it is used as the goal.  Column heads serve as joints of beams, and together with beams, they form a spatial framework, therefore it can be assumed that displacements are evenly distributed from their original place, i.e., they meet Laplace-Gaussian distribution. To bring this problem to numerical calculation, the column head could be simplified to one circle, which could be deconstructed into two parameters: diameter and a centroid (point).
To define an axis, at least two centroids are required, yet this axis is certainly determined by axiom. However, if there are three points or more, an axis is uncertainly determined by least-squares fitting. Consequently, it can be assumed that grids are rectangles, displacements of column heads are evenly distributed along each rebuild axis, and the variance of distances between current centroids and rebuild centroids of column heads has a minimum value. Therefore, Gauss sum instead of variance is expected to be minimum; thus, it is used as the goal.
The variance of a random variable X is the expected value of the squared deviation from the mean of X. That is, variance equals the average value of the Gauss sum. The squared root of the variance is the standard deviation, which is used as an evaluating indicator of the result.
Therefore, the problem is simplified as such: how can one generate all axes at the same time while the minimum value of Gauss sum is reached?
The methodology is shown as follows: 1. Use the Circle-FitPoints command to manually select related points from the point cloud of a column head and then automatically draw a circle of the section profile.
Repeat to go through all column heads ( Figure 3). 2. Launch Grasshopper, use the Curve component to collect all circles, use Project to obtain circles in the X-plane, use Area to obtain centroids of circles, use Deconstruct to obtain X, Y, and Z values, and use Sort List to sort centroids in a certain order ( Figure 4). 3. Use Number Slider to manually give a range for initial parameters including motions and rotations of axes, then connect them to the Genome of Galapagos as a candidate solution for the genetic algorithm ( Figure 5). The variance of a random variable X is the expected value of the squared deviation from the mean of X. That is, variance equals the average value of the Gauss sum. The squared root of the variance is the standard deviation, which is used as an evaluating indicator of the result.
Therefore, the problem is simplified as such: how can one generate all axes at the same time while the minimum value of Gauss sum is reached?
The methodology is shown as follows: 1. Use the Circle-FitPoints command to manually select related points from the point cloud of a column head and then automatically draw a circle of the section profile.
Repeat to go through all column heads ( Figure 3). 2. Launch Grasshopper, use the Curve component to collect all circles, use Project to obtain circles in the X-plane, use Area to obtain centroids of circles, use Deconstruct to obtain X, Y, and Z values, and use Sort List to sort centroids in a certain order ( Figure 4). 3. Use Number Slider to manually give a range for initial parameters including motions and rotations of axes, then connect them to the Genome of Galapagos as a candidate solution for the genetic algorithm ( Figure 5). 4. At this moment, a candidate solution and resulting axes are shown; the average distance is 144 mm ( Figure 6). 5. This solution could be optimized as long as a goal is given to the fitness of Galapagos, therefore Gauss sum is applied to calculate the applicability of axes, then the average distance between current centroids and rebuild centroids is applied to evaluate the result ( Figure 7). 6. Open the Galapagos, start the solver, and during the process, the result is displayed in real-time ( Figure 8). 7. When the result seems stable after 5 min of calculation, it could be output as the final solution, if the average distance is acceptable. Or wait until the minimum value is reached, though it may take longer (almost 1 h). The average distance is 68 mm ( Figure 9). result ( Figure 7). 6. Open the Galapagos, start the solver, and during the process, the result is displayed in real-time ( Figure 8). 7. When the result seems stable after 5 min of calculation, it could be output as the final solution, if the average distance is acceptable. Or wait until the minimum value is reached, though it may take longer (almost 1 h). The average distance is 68 mm (Figure 9).    As is shown in Figure 9, displacements of column heads including vectors (arrows) and lengths (numbers) are calculated, while all of the axes are generated. The solution process takes 3 h 16 min and stops automatically after reaching the minimum value.
This result shows that the manually rotated point cloud is very close to the orthogonal direction as the initial angle is 0.00002 degrees. The winding corridor has an average width of 3493 mm, and the four dimensions (3485 mm, 3488 mm, 3498 mm, 3500 mm) have a standard deviation of 6 mm. Therefore, the winding corridor could be regarded as equidistantly offset from the second circle of columns. Similarly, the whole grid system could be regarded as symmetric.
The components of the entire node-visible algorithm are shown in Figure 10.

Hypothesis and Verification
Ideally, if two tangent circles of 400 mm diameter arrayed in a vertical line are given, one axis should be generated; at this moment, each distance should be zero, the angle of the axis should be zero, and the aligned dimension between two horizontal axes should be also 400 mm. After the solution process, the result successfully verifies the hypothesis, although it takes 33 min to reach the goal ( Figure 11).    As is shown in Figure 9, displacements of column heads including vectors (arrows) and lengths (numbers) are calculated, while all of the axes are generated. The solution process takes 3 h 16 min and stops automatically after reaching the minimum value.
This result shows that the manually rotated point cloud is very close to the orthogonal direction as the initial angle is 0.00002 degrees. The winding corridor has an average width of 3493 mm, and the four dimensions (3485 mm, 3488 mm, 3498 mm, 3500 mm) have a standard deviation of 6 mm. Therefore, the winding corridor could be regarded as equidistantly offset from the second circle of columns. Similarly, the whole grid system  As is shown in Figure 9, displacements of column heads including vectors (arrows) and lengths (numbers) are calculated, while all of the axes are generated. The solution process takes 3 h 16 min and stops automatically after reaching the minimum value.
This result shows that the manually rotated point cloud is very close to the orthogonal direction as the initial angle is 0.00002 degrees. The winding corridor has an average width of 3493 mm, and the four dimensions (3485 mm, 3488 mm, 3498 mm, 3500 mm) have a standard deviation of 6 mm. Therefore, the winding corridor could be regarded as

Hypothesis and Verification
Ideally, if two tangent circles of 400 mm diameter arrayed in a vertical line are given, one axis should be generated; at this moment, each distance should be zero, the angle of the axis should be zero, and the aligned dimension between two horizontal axes should be also 400 mm. After the solution process, the result successfully verifies the hypothesis, although it takes 33 min to reach the goal (Figure 11). If there are four circles, equally distributed along the vertical axis, the result should restore the movement and rotation, and the average distance should be 200 mm. However, the result shows a wrong solution (Figure 12). The average distance is 182 mm, and the If there are four circles, equally distributed along the vertical axis, the result should restore the movement and rotation, and the average distance should be 200 mm. However, the result shows a wrong solution ( Figure 12). The average distance is 182 mm, and the initial angle is not zero. This is because four circles could form a rhombus, for which it is impossible to determine the only axis. If there are four circles, equally distributed along the vertical axis, the result should restore the movement and rotation, and the average distance should be 200 mm. However, the result shows a wrong solution ( Figure 12). The average distance is 182 mm, and the initial angle is not zero. This is because four circles could form a rhombus, for which it is impossible to determine the only axis. Therefore, the hypothesis is changed to eight circles, to form a symmetric grid system in which the horizontal dimension is 1600 mm and each vertical dimension is 800 mm, then these circles are distributed equally (200 mm) along each vertical axis.
The generated axes successfully restore the movement and rotation (Figure 13), although the solution process takes 54 min to reach the minimum value then stops automatically. However, when this grid system is asymmetric, the resulting axes have a smaller average distance (183 mm) than the initial pre-given axes (200 mm) and a smaller standard deviation (188 mm) than the initial pre-given axes (200 mm), but they fail to restore the original grids ( Figure 14). Therefore, the hypothesis is changed to eight circles, to form a symmetric grid system in which the horizontal dimension is 1600 mm and each vertical dimension is 800 mm, then these circles are distributed equally (200 mm) along each vertical axis.
The generated axes successfully restore the movement and rotation (Figure 13), although the solution process takes 54 min to reach the minimum value then stops automatically. However, when this grid system is asymmetric, the resulting axes have a smaller average distance (183 mm) than the initial pre-given axes (200 mm) and a smaller standard deviation (188 mm) than the initial pre-given axes (200 mm), but they fail to restore the original grids ( Figure 14).    Furthermore, both asymmetric and symmetric grid systems of six circles are tested. The only difference is that the top two circles are removed. Both resulting axes have a smaller average distance (177 mm) than the initial pre-given axes (290 mm and 202 mm), and both resulting axes have a smaller standard deviation (188 mm) than the initial pregiven axes (312 mm and 203 mm). The solver has reached the minimum value and stopped automatically (Figures 15 and 16). Furthermore, both asymmetric and symmetric grid systems of six circles are tested. The only difference is that the top two circles are removed. Both resulting axes have a smaller average distance (177 mm) than the initial pre-given axes (290 mm and 202 mm), and both resulting axes have a smaller standard deviation (188 mm) than the initial pregiven axes (312 mm and 203 mm). The solver has reached the minimum value and stopped automatically (Figures 15 and 16).  Additionally, to verify that an odd number of columns would have an inaccurate result, the hypothesis is changed to 10 circles. The generated axes successfully restore the rotation but fail to restore the movement. The solver has reached the minimum value and stopped automatically (Figures 17 and 18).  Additionally, to verify that an odd number of columns would have an inaccurate result, the hypothesis is changed to 10 circles. The generated axes successfully restore the rotation but fail to restore the movement. The solver has reached the minimum value and stopped automatically (Figures 17 and 18). Additionally, to verify that an odd number of columns would have an inaccurate result, the hypothesis is changed to 10 circles. The generated axes successfully restore the rotation but fail to restore the movement. The solver has reached the minimum value and stopped automatically (Figures 17 and 18).  This result shows that for linear-arrayed circles of a vertical axis, each column shares equal weight, so the three circles have more "attraction" than two circles, and the axis deviates to their side. The average distance is 192 mm, which is even more "ideal" than the original distribution of 200 mm, despite the axis not being the same as the original one. The standard deviation is 196 mm, which is closer to 200 mm, compared to 188 mm, of six columns. After the solution process is complete, the resulting axes of 10 columns are more accurate than six columns.

Results
The results of the hypothesis and verification are concluded in Table 1.  This result shows that for linear-arrayed circles of a vertical axis, each column shares equal weight, so the three circles have more "attraction" than two circles, and the axis deviates to their side. The average distance is 192 mm, which is even more "ideal" than the original distribution of 200 mm, despite the axis not being the same as the original one. The standard deviation is 196 mm, which is closer to 200 mm, compared to 188 mm, of six columns. After the solution process is complete, the resulting axes of 10 columns are more accurate than six columns.

Results
The results of the hypothesis and verification are concluded in Table 1. This result shows that for linear-arrayed circles of a vertical axis, each column shares equal weight, so the three circles have more "attraction" than two circles, and the axis deviates to their side. The average distance is 192 mm, which is even more "ideal" than the original distribution of 200 mm, despite the axis not being the same as the original one. The standard deviation is 196 mm, which is closer to 200 mm, compared to 188 mm, of six columns. After the solution process is complete, the resulting axes of 10 columns are more accurate than six columns.

Results
The results of the hypothesis and verification are concluded in Table 1. This method is proposed for the correction of an irregular grid system; the result could accurately restore the original axes under certain conditions. 1. At least two rows. 2. Each row has four symmetrically displaced columns or has more than four symmetrically displaced columns as long as the number is even (6, 8, 10, etc.).
Otherwise, the result could not accurately restore the original axes, but is still close to the original axis and still has potential possibilities for application, as it is repeatable and more convenient than manual calculation. The more columns a grid system has, the more accurate the result is.

Conclusions
This method could either output geometric models or numeric data. After the solution process, the resulting axes could be directly used for modeling. Alternatively, in a collaborative workflow including many sessions, to make data exchange convenient, only grid spacing and rotation angle are required. Both of them are numbers, and therefore easy to store in a general format such as .txt, .xlsx, etc. The grid spacing is used to create a grid system in BIM software such as Autodesk Revit and Bentley Openbuildings Designer. The rotation angle is used to rotate the point cloud to fit the grid system; thereafter, the point cloud is used as a reference during the modeling process.
Although this method aims to automate the solution process, regarding the diversity of architectures there are still manual operations. For example, removing some columns based on the standard column grids is a common practice, such as Kuiwen Pavilion. Not all axis intersections have columns; therefore, the corresponding column has to be removed from the algorithm.
Furthermore, this method could be applied to test whether these columns are designed to be inclined or not. First, the axes of column bases are generated. The centroid of rebuild axes on column bases could be recorded as C0, and the centroid of rebuild axes on column heads could be recorded as C1. The vector from C0 to C1 on the XY plane would indicate the movement between the column heads grid system and the column bases grid system. If the column heads grid system is moved to make the X and Y coordinates of C1 and C0 consistent, the result would indicate whether these columns are designed to be inclined or not (Figures 19 and 20).  If they are designed to be inclined, the outer eaves columns are inclined inward in the front and back eaves, with 10‰ of the column height; 8‰ of the column height in the two sides eaves, and the corner columns are inclined in both directions at the same time.
After the manual check, the result shows that the inclination of the columns almost conforms to a designed inclination, as the average distance is 50 mm, the column height is 5027 mm, and the ratio is almost 10‰. However, this result still needs further investigation to confirm, as all of the columns have inclinations, rather than just eaves columns. This practice has not been found in historical archives or architectural books. The cause of this special inclination needs further research to verify whether this construction method of column inclination is common in other Imperial-style architectures.    If they are designed to be inclined, the outer eaves columns are inclined inward in the front and back eaves, with 10‰ of the column height; 8‰ of the column height in the two sides eaves, and the corner columns are inclined in both directions at the same time.
After the manual check, the result shows that the inclination of the columns almost conforms to a designed inclination, as the average distance is 50 mm, the column height is 5027 mm, and the ratio is almost 10‰. However, this result still needs further investigation to confirm, as all of the columns have inclinations, rather than just eaves columns. This practice has not been found in historical archives or architectural books. The cause of this special inclination needs further research to verify whether this construction method of column inclination is common in other Imperial-style architectures.
Author Contributions: Conceptualization, X.W., C.W. and C.B.; Funding acquisition, C.W.; Investigation, X.W. and C.W.; Methodology, X.W., C.W. and C.B.; Project administration, C.W.; Resources, C.W.; Software, X.W.; Validation, X.W.; Visualization, X.W.; Writing-original draft preparation, X.W.; Writing-review and editing, X.W. and C.W. All authors have read and agreed to the published version of the manuscript.  If they are designed to be inclined, the outer eaves columns are inclined inward in the front and back eaves, with 10‰ of the column height; 8‰ of the column height in the two sides eaves, and the corner columns are inclined in both directions at the same time.
After the manual check, the result shows that the inclination of the columns almost conforms to a designed inclination, as the average distance is 50 mm, the column height is 5027 mm, and the ratio is almost 10‰. However, this result still needs further investigation to confirm, as all of the columns have inclinations, rather than just eaves columns. This practice has not been found in historical archives or architectural books. The cause of this special inclination needs further research to verify whether this construction method of column inclination is common in other Imperial-style architectures.

Data Availability Statement:
The source code of the node-visible algorithm is not openly available due to the consideration of data safety.