RooTri: A Simple and Robust Function to Approximate the Intersection Points of a 3D Scalar Field with an Arbitrarily Oriented Plane in MATLAB

: With the function RooTri() , we present a simple and robust calculation method for the approximation of the intersection points of a scalar ﬁeld given as an unstructured point cloud with a plane oriented arbitrarily in space. The point cloud is approximated to a surface consisting of triangles whose edges are used for computing the intersection points. The function contourc() of M ATLAB is taken as a reference. Our experiments show that the function contourc() produces outliers that deviate signiﬁcantly from the deﬁned nominal value, while the quality of the results produced by the function RooTri() increases with ﬁner resolution of the examined grid.


Introduction and Related Work
In the field of science and especially in engineering, the formulation of mathematical models for the treatment of concrete problems plays a key role.Often, however, the models are characterized by a high degree of complexity, which is why an analytical solution is only possible in rare cases.For this reason, one often has to rely on numerical methods to approximate the solution.A common application is the determination of roots, since these often represent characteristic points in the course of the function.For the numerical approximation of the roots of one-dimensional functions, the bisection method ( [1], p. 250), Newton's method ( [2], p. 243) or the method by R. P. Brent, who combines the secant method ( [1], p. 254) with the bisection method in [3], are prominent representatives.
On the other hand, the determination of roots of two-dimensional functions is much more complex.Here, the roots form a level set, which is represented as a contour line.In the context of our work, we focus on two-dimensional functions (scalar fields) of the form f : R 2 → R and assume that these are given as a point cloud P ∈ R n×3 , n ≥ 3.This has the advantage that even unstructured point clouds, which arise, for example, during the recording of measuring points, can be handled.In the following, we develop a simple and robust function for approximating the zeros of a point cloud with a plane oriented arbitrarily in space in MATLAB © .
The literature reveals that different methods are used to determine isolines or isosurfaces, such as the Marching Squares (2D) or Marching Cubes (3D) algorithm developed by W. E. Lorensen and H. E. Cline [4].It is mainly used in computer graphics and image processing.The core of the method is a divide-and-conquer approach, in which the data set at hand is divided into square grid cells and each cell is assigned a defined isovalue.The next step is to check whether the vertices are above or below the isovalue.By assigning a binary value, it is possible to determine which combination of corners of the cell corresponds to the isovalue, whereupon the polyline or polygon can be created.This is performed for all cells, from which the isoline or isosurface can be constructed.The method has been widely studied in research.Thus, C. Maple presents, in [5], an extension of the method for approximating the enclosed area or volume for use in room planning.Z. P. T. Sin and P. H. F. Hg use the Marching Cubes algorithm in [6] to generate planetary terrain in a spherical environment.In [7], L. Custodio et.al extend the Marching Cubes algorithm to ensure topological correctness.Another application can be found in [8].Here, X. Wang et al. develop an approach to implement the Marching Cubes algorithm based on edge growth.Closely related to the Marching Cubes algorithm is the Asymptotic Decider algorithm developed by G. M Nielson and B. Hamann, which generates isosurfaces from a given scalar field and provides topologically correct contours [9].A modified approach to robust and topologically correct triangulation is developed by R. Grosso in [10].In this context, A. M. Tushar and C. R. Johnson provide a framework based on the Asymptotic Decider algorithm for analyzing uncertainties related to isocontour generation [11].Other approaches regarding the determination of contour lines include the Octree-based algorithm [12] and the Span Space algorithm [13].
We note that the described algorithms basically aim to divide the considered scalar field into cells or cubes, to which a specific indicator is assigned in order to reconstruct the contour lines based on this information.These algorithms are primarily used in image processing, computer visualization, or medical imaging.In our research, we want to deviate from this approach and instead take an unstructured point cloud into account, using the information provided by the individual points directly.
With regard to commercially available software, MATLAB © offers contourc(), which can be used as a function to determine the contour lines of a two-dimensional function ( [14], p. 225).However, it is limited to a structured grid.For instance, O. Demirkaya et al. use it in terms of medical image processing ( [15], p. 86).In the further course of our investigations, we also use this function as a reference, since practice shows that it is often applied.

Modeling and Implementation
The following section covers, in its first part, the description and mathematical modeling of the underlying problem.Based on this, we then develop the function and implement it in MATLAB © .

Mathematical Modeling
In the first step, we wish to generate a surface from a given point cloud P ∈ R n×3 containing elements p i ∈ P, i ∈ I, with I ⊂ N being the associated index set.An element p i further consists of its specific components and is defined as From the definition of the scalar field described in Section 1, it follows that each pair of values [x i , y i ] is assigned exactly one real number z i .This means that we are able to remove the z i component first and project the remaining component pairs [x i , y i ] onto the xy-plane, which is qualitatively illustrated in Figure 1.From this procedure, we obtain a subset P ∈ R n×2 containing the two-dimensional points denoted by p i = [x i , y i ], i ∈ I.These can then be utilized to perform DELAUNAY triangulation to produce triangles that do not intersect.Moreover, the indices of points that form a triangle are known, and this information will be of importance in the further course of the work.Note that DELAUNAY triangulation is a common method for generating meshes in finite element analysis ( [16], p. 394).A detailed description of the method can be found in ( [17], p. 199).The vertices of the triangles are then reassigned to the corresponding z i components, and we obtain the surface approximated using triangulation.This approach offers the advantage that the basic structure of the triangles is preserved, and we additionally obtain a set K ∈ R m×3 containing all triangles k defined by the indices of the points.In the second step, the intersection points of the surface with an arbitrarily oriented plane are to be calculated.For this purpose, we first consider a random triangle k ∈ K of the approximated surface.The idea now is to use vertices v k j , j ∈ {1, 2, 3} of the triangle for generating vectors t k j , j ∈ {1, 2, 3}, which are examined to determine whether they intersect the plane.With λ ∈ [0, 1], these are to be computed as follows: Consider the plane given in its coordinate form On this basis, intersection point s k j of a vector t k j with this plane is to be determined first.For this purpose, its components are inserted into the plane equation defined by Equation ( 1) and solving for λ yields Then, the λ given by Equation ( 2) is inserted back into the equation for t k j to compute intersection point s k j , and it is subsequently to be verified whether the latter is an element of the interval defined by the respective vector.Therefore, we set the system of equations and solve it for the new λ * .In case s k j is an element of the specific interval and is thus a relevant intersection point, λ * must be in turn an element of the interval [0, 1].From this follows the sufficient condition It is easy to see that the number of possible intersection points depends on the number of triangles and thus on the resolution of the approximated surface.Therefore, to make the evaluation more efficient, we propose to exploit the information provided by the triangle and define additional auxiliary points a k j , j ∈ {1, 2, 3} on the midpoints of the edges.Since the vectors are already defined, these points can be easily computed by setting λ = 0.5.With the help of the auxiliary points, new vectors can now be generated and are used to determine further intersection points according to the same scheme.Hence, the number of intersection points increases since for each triangle; in total, nine vectors are evaluated, while the resolution of the surface remains unchanged.In this context, Figure 2 qualitatively shows our approach for a random triangle k ∈ K intersecting a plane.In addition, Figure 3 illustrates common types of offset planes with the corresponding values for the plane parameters.

Note on the Determination of Plane Parameters
We present a brief tutorial on how to determine the plane parameters and consider three appropriate points x i ∈ R 3 , i ∈ {1, 2, 3} that can be used to construct the plane.First, we define the parameter form with x 1 as the support vector and obtain In this context, we further require that the direction vectors are not collinear, which means that the condition must be satisfied to ensure that they are not multiples of each other and thus do not point in the same direction.Otherwise, it would not be possible to span a plane.Then, we determine normal vector n using the cross product, i.e., and insert it into the coordinate form in Equation ( 1).The remaining parameter, d, is obtained by inserting support vector x 1 .Finally, the parameters result in This approach can be used to determine the parameters of an inclined plane.

Numerical Implementation and Code Description
The method described in the previous section is now to be implemented as a function called RooTri(arg1, arg2, arg3, arg4, arg5) in the MATLAB © software environment.A total of five arguments are provided as input variables to the function and are listed in Table 1.In the first step, the input variables are read to define the plane and extract the x and y components of the elements of P. After that, DELAUNAY triangulation is performed using the function delaunay() provided in MATLAB © .This produces a matrix K ∈ R m×3 containing the indices of those points from P that form a triangle.Subsequently, the vectors and intersection points are determined for all elements of the matrix, which are also checked to verify whether they are elements of the interval defined by the vector in question using Equation (4).If this is true, these intersection points are stored in an output matrix P ∈ R p×3 , which contains the components of all valid intersection points.This process is repeated until all triangles of matrix K have been checked.In case no intersection point is found, the matrix contains no element, i.e., P = ∅.Here, Figure 4 summarizes the general procedure of the algorithm.In the course of increasing the performance of the algorithm, we refrain from using for loops in the implementation but instead build the code vectorized.This has the advantage that the code can be executed faster, which is especially important if the function has to be called frequently.However, the general principle of operation based on Figure 4 remains the same.The entire code can be found in Appendix A.

Experiments and Validation
Subsequent to the implementation, the function RooTri() shall be compared with the function contourc() provided in MATLAB © , which is usually employed for computation of low-level contour matrices and can be used to determine the roots of a two-dimensional function.Due to the abundance of comparison possibilities, it is obvious that the determination of statistical parameters requires a methodologically well-founded procedure.Furthermore, our experiments were performed under the system specifications listed in Table 2 below.

Design of Experiments
We choose the method of the so-called Latin Hypercube Design (LHD), a wellestablished procedure in the field of Design of Experiments (DoE).Using the Latin Hypercube Design, a test field can be created in such a way that the variance of the global mean of the variables under consideration is minimized.Compared with other methods, such as the Monte Carlo method, it provides better results in terms of minimizing variance with the same number of test points.Thereby, a q × r matrix called L q×r , whose columns consist of an arbitrary permutation of the numbers {1, 2, 3, . . ., q}, is formed.In the framework of the so-called Latin Hypercube Sampling (LHS), a random number from [0, 1) is subtracted from each value of the LHD and then divided by q, thus normalizing the entire test field ( [18], p. 155).In our investigation, we use the function lhsdesign() from MATLAB © .Furthermore, it is necessary to limit the number of possible combinations and thus the dimensions of the test space to a feasible level.In this connection, J. Loeppky et al. propose, in [19], to estimate the number of experiments (n t ) considering a given number of dimensions n D using the following formula: Usually, the variables considered in the LHD are continuously distributed.Nevertheless, we also want to consider different types of signals and thus include several functions as variables that are added randomly.Hence, to define the parameters to be varied, we first specify four functions in total that intersect a plane with offset d parallel to the xy-plane; see Table 3.Here, PARABOLA is a common two-dimensional function and is depicted in Figure 5. On the other hand, the PEAKS function depicted in Figure 6 is a standard function of MATLAB © that is characterized by prominent minima and maxima.Furthermore, we investigate the ROSENBROCK function [20] as well as HIMMELBLAU's function [21] (see Figures 7 and 8).Note that these functions are generally applied as test functions for evaluating algorithms for solving nonlinear optimization problems.As the next parameter to be varied, we choose variable d, which determines the offset of the plane.We want to underline at this point that the RooTri() function is also able to consider arbitrarily oriented planes; however, we limit the analysis to parameter d, since the contourc() function can only determine intersections with parallel planes.Thus, we are able to compare the two functions with each other.The last parameter refers to length 2l of the interval on which the functions are to be examined.We define this in such a way that it is divided into a certain number of equidistant increments and extends over a range from −l to l in both the xand y-directions using the function linspace().This gives us a total of three parameters that are used to generate the test field while we vary (half) length l of the interval in a range from 1 to 10 and offset parameter d in a range from −5 to 5. According to Equation (5), 30 experiments must, therefore, be carried out and are summarized in Table 4.In this context, we additionally introduce another quantity, denoted by α, that relates the length of the interval (2l) to the number of increments (n α ).
Taking α and 2l as known, Equation ( 6) can be used to compute the required increment number (n α ), which ensures that the experiments are performed under the same conditions with respect to the resolution of the interval.From this, we initially obtain a regular grid on which the functions are to be evaluated.However, it is known from practice that, for instance, when recording measuring points, these values to obtain a regular grid cannot be set exactly.Instead, an error is produced, which leads to a distortion of the regular grid and results in unstructured data points.We additionally consider such an effect and assume a percentage deviation ε from the ideal value (ε = 0 just corresponds to the regular grid).This describes a type of tolerance band in which the values can deviate.With the increase in ε, the grid is distorted, which is qualitatively shown in Figure 9.

Name Equation
PARABOLA function In the course of numerical experiments, we mainly want to focus on the following aspects for both the contourc() and RooTri() functions: The quality ∆ f i (α, ε), i ∈ {1, . . ., 4} refers to the function value that is obtained by reinserting the calculated roots into the respective test function.The smaller the deviation from d is, the more accurately the algorithm solves the problem.Considering a computed intersection point s, the above is to be calculated as follows: Furthermore, we want to evaluate the number of intersection points, n r , computed by the RooTri() function and n c of the contourc() function.As a rule, the contour lines can be constructed more effectively with a greater number of intersection points.This is especially the situation as soon as the isolines are close together.In case of large data volumes, both required computing time t r , t c and required memory q r , q c are of significant importance.For this reason, these parameters are also to be evaluated.

Evaluation of Computed Results
As part of the evaluation, all experiments from Table 4 are performed for a given combination (α, ε).Here, α is varied on an interval between 0.025 and 0.250, and ε, on an interval between 0.000 and 0.250, each divided into 25 equidistant increments, resulting in a total of 30 × 625 tests.In turn, the results are used to determine the arithmetic mean and the empirical standard deviation of the corresponding variables under consideration.

Quality of Computed Intersection Points
The evaluation of the quality of the calculated results for both functions is shown in Figure 10 in logarithmic scale.In this connection, the arithmetic mean, µ, as well as the empirical standard deviation, σ, were calculated by applying the equations Based on Figure 10, we note that with finer resolution of the grid (i.e., with smaller α), the arithmetic mean of the percentage deviation from the nominal value decreases when using the RooTri() function.The same applies to the evolution of the standard deviation.For example, taking into account a regular grid (ε = 0) with α = 0.025, we obtain an average percent deviation from the nominal value of 2.49% ± 3.46%.It can be seen that the percentage deviation converges to zero with smaller α.This basic tendency remains the same with the increase in ε.
In contrast, we observe an opposite course when applying the contourc() function.Paradoxically, the percentage deviation from the nominal value increases with finer resolution of the grid.Investigations show that the function produces outliers that deviate from the nominal value to a high degree.Thus, the standard deviation also increases with the decrease in α.Again, the tendency remains the same with the variation in ε.To illustrate this fact, we demonstrate the performance of the contourc() function by first defining a further test function: We evaluate f 5 (x, y) on the interval [−25, 25] × [−25, 25] for different α and determine the intersection points at which f 5 (x, y) = 0 applies; see Figure 11.It can be clearly seen that here, the contouc() function also produces outliers that deviate significantly from the nominal value.For example, the plotted outlier for α = 2/3 at [0.0, 83.0] leads to a function value of approximately f 5 (0.0, 83.0) ≈ 0.48, while the observed effect increases with smaller α.An exemplary application can be found in Appendix B. At this point, it should be emphasized that although these outliers can be made visible in the course of a graphical representation, they can only be identified with great effort when simply looking at the generated results in the form of an output matrix.This leads to the necessity of an additional processing of the data.

Number of Computed Roots
Figure 12 shows the arithmetic mean of the ratio of the calculated intersection points, n r /n c .From the diagram, it can be concluded that for a regular grid, the ratio is almost constant, and on average, six times more intersection points are calculated using the RooTri() function.However, with the increase in error ε, the ratio decreases.From this, we deduce that the number of calculated intersection points depends on the quality of the grid.

Time Consumption
The time for computation is determined by setting a timestamp before and after the function call is completed in each case.Thereby, we evaluate all performed experiments and obtain a ratio of t c /t r = 1.1285 ± 0.087.From this, we find that the contourc() function takes 12.85% more time on average.This is mainly due to the fact that when considering a distorted grid, the griddata() function must also be used to interpolate the distorted grid back to a regular grid.

Estimation of Required Memory
The exact determination of the memory requirement is complicated; for this reason, it is to be estimated over the elements to be processed.Here, we again assume that a point cloud P ∈ R n×3 is given.For the application of RooTri(), therefore, 3n elements must be processed.In the case of contourc(), in addition to the 3n elements, 3 √ n √ n elements must be considered due to the matrices required for the griddata() function.Thus, we obtain the ratio to estimate the required memory.On this basis, we can formulate the statement that the contourc() function takes up twice the memory used by the Rootri() function for the same task.

Limits of Application
The quality of the results depends, in particular, on the degree of nonlinearity of the underlying function.In Figure 13, this circumstance is qualitatively illustrated for a true function with two points p 1 and p 2 from a given point cloud.In case of a moderate change rate of the true function (Figure 13a), the approximation error to the true intersection remains small, but if the behavior of the function between these support points is strongly nonlinear, the deviations may also increase; see Figure 13b.In such a situation, it is recommended to increase the resolution of the considered grid.Based on the evaluation, we can draw the conclusion that the RooTri() function is significantly more robust than the contourc() function in terms of the quality of the computed results.Furthermore, the calculation is performed faster, and less memory is required, since fewer elements have to be processed.At this point, it should be further emphasized that the quality of the roots could only be determined because the test functions were known in advance.If these are not known, there is no possibility to detect outliers, since the control possibility is omitted.Rather, one is then dependent on the manual preparation of the data to identify the outliers.

Practical Application
Following the evaluation, the application of our function will be demonstrated with practical examples to highlight the wide range of possible fields of usage.

Exemplary Function Operation
As a generic example for the application of the function RooTri(), we consider a point cloud P that is generated from the function and is available as a file called RooTriExample.txt.It is first defined as P and loaded using the readmatrix() function.In the following step, the parameters describing the plane are defined and passed to our function alongside point cloud P, i.e., we execute RooTri(P,a,b,c,d).Here, Figure 14 shows, on the one hand, point cloud P and, on the other hand, the resulting intersection points at different parameters for the description of the respective intersection plane that we arbitrarily selected for our example (regarding a certain definition of the intersection plane, see Section 2).It is clearly visible that the function RooTri() can also be used for arbitrarily oriented planes to determine the intersection points.The corresponding example script can be found in Appendix C.

Measuring Electrical Machines
In our next example, we consider the measurement of a permanently excited synchronous machine for which real measured values were recorded and appropriately normalized.A qualitative illustration of the machine is depicted in Figure 15.Torque M of the machine depends on current i, which is composed of components i d and i q , so the relation applies.Component i d of the current generates a magnetic field in the direction of the flux linkage of the permanent magnet, and the i q component generates a magnetic field that is rotated by 90 • with respect to the flux linkage [22].Within the scope of the investigation, a characteristic diagram of the machine that contains the equipotential lines of constant torques is to be created on the basis of the measured values.1) describes, in this case, the value of the respective torque M, which is varied in a range from −0.9 max{|M|} to 0.9 max{|M|}.The function is thus executed by calling RooTri(P,0,0,1,d) for each ratio M/ max{|M|} as offset parameter d.In Table 5, the results with respect to the number of calculated zeros (n) and the required computation time (t) are summarized.Using the calculated combinations [i d / max{|i d |}, i q / max{|i q |}] to obtain the respective torque ratio, the characteristic diagram can then be generated; see Figure 16.Note that due to the number of data, not all results were plotted; however, the qualitative course of the equipotential lines is clearly recognizable.The function offers the further advantage that the points of the isolines are explicitly available.These working points can then be approached directly for later operation.

Structural Optimization
As a further use case of the RooTri() function, we consider a typical problem from the area of so-called topology optimization, a subarea of structural optimization.The focus is on the search for a suitable material distribution in a given design domain with defined boundary conditions and external loads in order to fulfill a certain objective criterion.As a benchmark problem, compliance is often minimized.Figure 17a shows such a problem.Here, a linear elastic body that has corresponding Young's modulus E, density ρ and POISSON's ratio ν is considered.On the left side, the body is clamped, and on the right side, a force F acts.The extensions in the xand y-directions are marked with 2L and 2H, respectively.To determine the material distribution, different approaches are investigated in research, and a famous representative is the so-called level set method.Here, it is assumed that the component boundaries of the structure can be interpreted as the zeros of a higherdimensional function, the φ-function.At the component boundary, then, the condition φ(x, y) = 0 holds.For example, Z. Zhuang et al. apply the level set method in [23] and also develop a so-called adaptive mesh.This method can be used to discretize the component boundaries more finely.However, the complexity of the optimization problem requires to develop the φ-function iteratively to the local optimum in the framework.Due to the numerics, the φ-function is then available as a point cloud; in Figure 17b, this is shown for a certain iteration step.The task, now, is to intersect the point cloud with the xy-plane in order to determine the optimized contour.We, therefore, pass point cloud P depicted in Figure 17b to the RooTri() function and set parameter d to 0, i.e., we call RooTri(P,0,0,1,0).As a result, we finally obtain the associated component contour in the respective iteration step; see Figure 18.At this point, it should be noted that it is particularly useful, in the case of topology optimization, to obtain as many zeros as possible.Thus, the component contour can be described more precisely, which in turn is important for the discretization of the component in the context of finite element analysis.

Further Impact
The practical examples shown above illustrate two application options from the engineering field; however, the use of the RooTri() function is not limited to these examples, as it can also be used in many different ways.A classic field of application for determining isolines can be found in geography.Here, land surfaces are surveyed in order to extract elevation levels from geographical information.Scalar fields are also the subject of research in meteorology.Atmospheric variables such as pressure and temperature distributions or humidity are examined as parameters that can be used to create weather maps.Other areas of application can be found in biology or chemistry when examining substance concentrations, for instance.

Conclusions and Outlook
In this paper, we present a simple and robust function named RooTri() for determining the roots of a two-dimensional scalar field with a plane arbitrarily oriented in space.One major advantage is that our function is able to handle even unstructured data in comparison to the contourc() function provided in MATLAB © , which is used as reference in this work.To evaluate the function, a test field is generated based on proven statistical methods, and a total of four different test functions are examined.Our experiments show that with finer resolution of the interval to be investigated, the percentage deviation from the nominal value decreases and converges to zero.On the other hand, the percentage deviation due to outliers produced by the function contourc() even increases with finer resolution of the grid.Based on our results, we further find that compared with the contourc() function, our function uses half as much memory and works 12.85% faster on average.In addition, for an ideal grid, approximately six times more intersection points are calculated.Finally, the application of our function is demonstrated using practical examples from the field of electrical machines and structural optimization.As a result of the findings presented, we will use the RooTri() function for our further research work and have already been able to implement it successfully in more complex algorithms (for example, for automated measurement data evaluation or topology optimization) as a function module.In particular, the robust handling of unstructured data and the speed of the evaluation are decisive factors.
Concluding, we would like to highlight that the RooTri() function can be used for a variety of problems in which scalar fields in the form of point clouds play a role.Accordingly, it is also a cross-domain tool and is not exclusively limited to engineering problems.Therefore, we hope that researchers from a wide range of fields will use our function for their work.For future activities, the code should still be transferred to other programming languages (such as PYTHON) and made available to the community.% p l o t r e s u l t s s c a t t e r 3 ( ipmat ( : , 1 ) , ipmat ( : , 2 ) , ipmat ( : , 3 ) , 1 , ' k ' ) a x i s equal t i t l e ( ' I n t e r s e c t i o n p o i n t s ' ) x l a b e l ( ' x ' ) ; y l a b e l ( ' y ' ) ; z l a b e l ( ' z ' )

Figure 1 .
Figure 1.DELAUNAY triangulation and surface approximation; given point cloud P and predefined plane, projection of points to origin plane to obtain P , triangulation to obtain set of triangles K, adding z component to vertices of each triangle k ∈ K.

Figure 2 .
Figure 2. Determination of intersection points for a given triangle k ∈ K.

Figure 3 .
Figure 3. Origin planes with constant offset and corresponding plane parameters.

Figure 9 .
Figure 9. Regular grid (a) and distorted grid (b) due to incorrect settings.

Figure 10 .
Figure 10.Arithmetic mean and according standard deviation of the quality.

Figure 12 .
Figure 12.Arithmetic mean of the ratio of computed roots, n r /n c .

Figure 13 .
Figure 13.Influence of nonlinearities; (a) moderate nonlinear course of the true function, (b) strongly nonlinear course of the true function.Due to numerical inaccuracies, it may be further advisable in case of parallel plane sections to check the results again afterwards to see whether a pair of values exist twice in the section plane.For example, for planes parallel to the xy-plane with offset d, two pairs [x * , y * , d] and [x * , y * , d ± η] with η 1 can occur.One option is to only print the first two values, i.e., x and y components, since the corresponding value for d is known, and use the unique() command to remove duplicate pairs to avoid ambiguities.Based on the evaluation, we can draw the conclusion that the RooTri() function is significantly more robust than the contourc() function in terms of the quality of the computed results.Furthermore, the calculation is performed faster, and less memory is

Figure 15 .
Figure 15.Qualitative illustration of a permanently excited synchronous machine.A total of 93,853 measured values p in the form p = [i d / max{|i d |}, i q / max{|i q |}, M/ max{|M|}] (14) are available as unstructured data in a point cloud P and are passed to the RooTri() function.Parameter d of the plane equation given in Equation (1) describes, in this case, the value of the respective torque M, which is varied in a range from −0.9 max{|M|} to 0.9 max{|M|}.The function is thus executed by calling RooTri(P,0,0,1,d) for each ratio M/ max{|M|} as offset parameter d.In Table5, the results with respect to the number of calculated zeros (n) and the required computation time (t) are summarized.Using the calculated combinations [i d / max{|i d |}, i q / max{|i q |}] to obtain the respective torque ratio, the characteristic diagram can then be generated; see Figure16.Note that due to the number of data, not all results were plotted; however, the qualitative course of the equipotential lines is clearly recognizable.The function offers the further advantage that the points of the isolines are explicitly available.These working points can then be approached directly for later operation.

Figure 17 .
Figure 17.Benchmark problem in topology optimization; (a) problem definition, (b) numerically determined values of the higher-dimensional φ-function in a certain iteration step.

Figure 18 .
Figure 18.Resulting structure in a certain iteration step; component boundary represented by intersection points.

Table 3 .
Mathematical formulation of the considered test functions.

Table 4 .
Test field for numerical experiments.
Figure 16.Determination of equipotential lines of constant torques.