Development of the Object-Oriented Dynamic Simulation Models Using Visual C + + Freeware

The paper mostly focuses on the methodological and programming aspects of developing a versatile desktop framework to provide the available basis for the high-performance simulation of dynamical models of different kinds and for diverse applications. So the paper gives some basic structure for creating a dynamical simulation model in C++ which is built on the Win32 platform with an interactive multiwindow interface and uses the lightweight Visual C++ Express as a free integrated development environment.The resultant simulation framework could be a more acceptable alternative to other solutions developed on the basis of commercial tools like Borland C++ or Visual C++ Professional, not to mention the domain specific languages and more specialized readymade software such as Matlab, Simulink, and Modelica. This approach seems to be justified in the case of complex research objectoriented dynamical models having nonstandard structure, relationships, algorithms, and solvers, as it allows developing solutions of high flexibility. The essence of the model framework is shown using a case study of simulation of moving charged particles in the electrostatic field. The simulation model possesses the necessary visualization and control features such as an interactive input, real time graphical and text output, start, stop, and rate control.


Introduction
This paper is primarily intended for the researchers who focus on the development of scientific continuous and hybrid simulation models of dynamical systems and have programming skills in C++ and Windows.The dynamical systems evolving with the flow of time are described by sets of differential equations, ordinary or partial ones.Such systems are the main part of diverse modelling applications belonging to different scientific and engineering domains as high-energy and accelerator physics (see, e.g., Hajari et al. [1]), ship navigation (Sandaruwan et al. [2]), electric power systems (Yusop et al. [3]), mechatronics and robotics (Ferretti et al. [4]), aerospace (Kozynchenko [5]), process control (Sivakumaran and Radhakrishnan [6]), and so forth.These applications, being characterized by advanced theoretical level of underlying mathematical models, can be multidomain and have rather complex structures including frequent decision-making and intelligent control units.
The problem of choice of programming language and simulation environment which could be most suitable for developing such complex and unconventional models from the viewpoint of their universality, flexibility, effectiveness, and interactive features comes to the fore.There exists specialized ready-made software, such as Maple, Mathematica, MathCAD, and above-mentioned Matlab and Simulink, which provides a range of opportunities sufficient for a number of applications of average complexity, especially for the models consisting of standard subcomponents.On the other hand, a number of domain-specific languages have been developed for the simulation purposes, such as the one described by Ewald and Uhrmacher [7].However, it seems rather difficult to find the adequate tools for more complicated models.Even the use of the declarative languages like Modelica [8] is also questionable in this case, because the solution based on imperative (and compiled) languages, the C++ first of all, usually provides better effectiveness and flexibility.An example of applying the imperative language Object Pascal (Delphi) to the problem of simulation and optimization of the particle beam emitter has been shown by Kozynchenko and Svistunov [9].However, the rapid development environments like Delphi or its counterpart C++ Builder are proprietary software, whereas the use of open or free tools would be more suitable in case of lack of funding.
It should be noted that the problem of numerical simulation and visualization using Microsoft Visual C++.NET (which is a proprietary software) has been considered by Salleh et al. [10,11].The main feature of the approach discussed in the book is the support of MFC (Microsoft Foundation Classes) library in developing the graphical user interface applications.However, the essential issue of simulation of dynamical systems with the real time graphical output and simulation rate control is not addressed in the book.
To sum up, we can say that the problem of developing an appropriate flexible interactive graphical interface and framework for complex object-oriented dynamical simulation models still remains actual.This interface and framework should provide the following desirable features and properties for the dynamical simulation models focused on complex scientific and engineering problems: (1) Usage of efficient (presumably compiled), flexible, and widespread programming language, as well as a mainstream integrated development environment (IDE).
(2) Ability to draw graphs in progress, tracing out the change of state variables and other parameters with respect to the simulation time.
(3) Ability to interrupt the program execution, to resume and continue it, and to change the simulation rate over the time.
(4) Ability to make required inputs in the course of simulation, using several windows if necessary.
(5) Ability to display the selected numerical output in progress and simultaneously with displaying the graphical output.( 6) Usage of open-source and/or nonproprietary software.
In the paper, we try to show an example of a simulation framework that possesses the above features.The modelling application is taken from electrodynamics domain and describes the motion of charged particles in the field of an electrostatic lens.The model is developed in the standard unmanaged C++ using the Visual C++ Express 2010 Edition as a development environment.When creating a new project, the "Empty Project" template has been selected for simplicity, because it only fixes the file structure and, for the rest, allows the developer to create the project of own type.One can choose the "Win32 Project" template as well, and on the "Application settings" page the option "Empty project" should be selected.We focus on a Win32-based application, that is, an application that handles Windows messages directly using message loops.The pure Win32 C++ code is used to create the project from scratch.Such an alternative as writing a code with the help of the MFC library is not available, because the MFC is not supported in the Express Edition.The designed program relates to so-called dialogbased applications which use a dialog box as a main window.The project consists of the set of header files, field.h,particle.h,rk4 4eq.h, and resource.h,two resource files, script.rcand script2.rc,and a source file, main.cpp.

Structure of the Simulation Model
First of all, the general outline of the simulation model to be studied is provided below.The diagram shown in Figure 1 displays the class and function relationships and is helpful in introducing the subject.The model has an ordinary structure of the desktop dialog-based Win32 two-threaded programs and includes the following global functions and classes: (i) The function WinMain as an entry point at the main thread.
(ii) The application-defined callback functions ModalDl-gProc and ModelessDlgProc that process the messages in the modal and modeless dialog boxes.
(iii) The application-defined callback function ThreadProc that is executed by the second thread and performs the application-specific algorithms.
(iv) The application classes Field and Particle describing the motion of a charged particle in the field of the electrostatic immersion lens.
(v) The function template RK4 4eq which defines a Runge-Kutta solver carrying out one integration step for a set of ordinary differential equations (ODEs).
(  fragments do not contain some features being important but not affecting the program logic, like exception handling, private class sections, and set and get functions.

Description of the Modeling Application
Below we consider in more detail a physical model to be simulated.This is a two-dimensional model that includes a number of hydrogen nuclei (protons) moving in the axialsymmetric electrostatic field of the immersion lens (see, e.g., Humphries Jr. [12], §6.6).The model consists of two classes, Field and Particle, whose objects are created in the program.the diameter, of the immersion lens.These potentials are calculated by solving the Laplace equation with the use of the finite difference method and Liebmann's iteration process (see, e.g., Veerarajan and Ramachandran [13], Chapter 11).
The class Field contains also arrays of pairs (  ,   ) of the coordinates of the th point belonging to one of the equipotential lines to be plotted for potentials −10, −25, . . ., −990 V.A data member zMax being equal to the length of the immersion lens sets the maximal distance of particle motion along axis.This variable is used in the second thread in the function ThreadProc as a termination condition of simulation for the given particle (see Algorithm 9, line 21.9).

Class Particle. Another application-related class is the class
Particle defined in the header file particle.h(Algorithm 2).The class contains all variables and constants being necessary for calculating the motion of a proton in the electrostatic field.These are time , coordinates , , velocity V and its components V  and V  , kinetic energy , proton mass  and charge  presented as static constants, components of electrostatic force   ,   , the dimensions of the lens (lengthOfElectrode and lengthOfGap), a set of variables (such as uLeft and yCellTop) used in calculating the forces in the helper function FzFy, a field grid spacing  = , and a control variable kVy which is used in the Particle's constructor to set the initial values of velocity components V  and V  under the given energy  0 by the formulas: The Boolean flag variable stop initialized as false shows whether the particle moves inside the lens; it becomes true when it reaches the rightmost edge of the lens at   ≥  and so it falls out of the further simulation.Additionally, the class includes a nested structure Data, objects of which keep the information about the current particle's state and are placed in the STL vector container history for holding the particle's state during the simulation.The data members are initialized in the combined (i.e., both default and parameterized) constructor.Other member functions are four righthand sides of the ODE's set: DzDt, DyDt, DvzDt, and DvyDt.The class also contains a helper function FzFy computing the values of electrostatic forces   and   used in the member functions DvzDt and DvyDt.

Description of the Two-Threaded Simulation Engine
Algorithm 3 lists the source file main.cpp that contains the function WinMain, as well as a set of global variables, constants, and function declarations (lines 10.3-27.3).Among them there are handles (descriptors) to dialog boxes hModalDlg and hModelessDlg, a handle hThread to a thread which carries out the model computations, control variables pause and sleep responsible for interrupting, resuming, and changing the rate of the simulation, constants height and width that determine the size of a plot rectangle where the particle trajectories, electrodes, and the field of the lens will be drawn, and scale factors my, mz, and kzy used in displaying the field.Also noteworthy is an STL vector container vpp which is used for keeping a set of Particle's objects taking part in the simulation.A global variable particleModelessBox is introduced to provide the particle's data numerical output at the modeless dialog window, along with a global variable particleListBox that indicates the selected particle whose numerical data is expected to be displayed in the list box placed at the modal dialog window.The function WinMain first creates a dynamically allocated object of the class Field, initializing the static pointer Particle::pF (line 31.3).Simultaneously, the electrostatic field is computed in the invoked constructor of the class Field.Furthermore, WinMain creates a modal dialog box from the dialog box template resource, using the function DialogBox.The corresponding resource-definition script (file script.rc) is listed in Algorithm 4. It contains most of resources relating to the modal dialog box and forming the graphical user interface (GUI), such as Windows controls (child windows): buttons, edit controls, text controls, a list box, and a trackbar.In Figures 2(a) and 2(b), a screenshot of the dialog boxes in the process of simulation is shown (a) together with a part of the modal dialog box with the most important GUI controls and associated resources' IDs (b).The function WinMain includes also a message loop that retrieves the Windows messages from the message queue (lines 33.3-38.3).If the messages are intended for the modeless dialog box, they should not be processed at the program's message loop.In this case, the function IsDialogMessage sends the messages directly to the modeless dialog box avoiding the program's message loop, as is done using the if statement in line 34.3.The function's body begins with definitions of the variables used in designing the graphical interface.These are a handle hDC to the device context, a structure defining a set of graphic objects, with some of them also defined here (pens, brushes, and fonts), a variable psPaint of type PAINTSTRUCT keeping the information necessary for drawing in a client area of the modal dialog box, a variable rect of type RECT used to keep coordinates of the client area, and temporary text buffers.It should be noted that most of the variables above must be declared as static ones; otherwise the graphical information will not be displayed properly.The basic part of the dialog box procedure is a switch statement with the controlling expression.The function is intended to be processing the Windows messages sent to the dialog box and passed through the second formal parameter message.The first message is WM INITDIALOG which is sent by the system immediately before displaying the dialog box and carries out the required initializations.In the case of sending this message, the following actions are made:

Description of the
(i) Creation of the second thread with the handle hThread in the suspended state (lines 27.5-32.5).
(ii) Initialization of the global handle to the modal dialog box hModalDlg (line 33.5) in order to be used further in the thread procedure as a control variable.
(iii) Creation of the modeless dialog box with the global handle hModelessDlg, which can display more detailed console-like output for the particles selected during the simulation (lines 34.5-35.5).
(iv) Completion of creating an edit box (IDC EDIT ParticleModelessBox) with an up-down (spin) control (IDC UpDown) which is sometimes referred to as a "spinner" control; this is done by calling the functions SendMessage to set the edit box as a buddy one and to set the range and the initial zero value (lines 36.5-40.5).The size of this container, that is, the total number  of particles involved in the simulation, is then passed to the edit window with an identifier IDS EDIT PartNum to be displayed (lines 13.6-14.6).

Modelling and Simulation in Engineering
The considered simulation model allows displaying not only graphs but also numerical data in the console-like fashion.Therefore, after at least one particle is created, the user may choose a particle whose state data could be displayed both at the above-mentioned modeless dialog box with the handle hModelessDlg and in the list box placed at the modal dialog box.It is accomplished in the edit box with the up-down control, where the required particle is chosen by the spin button.When the user clicks the button "Confirm," the numeric input, that is, the particle's number, is retrieved by the function GetDlgItemInt and is assigned to the variable particleModelessBox (lines 18.6-19.6).The particle chosen for the first time, that is, under the first click of the button "Confirm," is selected to be displayed in the list box.The number of this particle kept in the variable particleModelessBox is assigned to the variable particleListBox in line 23.6 (and made visible in line 24.6) which is then used in both threads, remaining unchanged during the simulation.It should be noted that the user can change the particle, whose data are displayed at the modeless dialog box, at any time, using the above procedure.
The modal dialog box contains also two buttons "Run" and "Pause" placing an important part in controlling the simulation process.When the user clicks the button "Run" the program checks if the global variable pause equals true.If so, the function ResumeThread is called and the value of pause is reversed to false (lines 27.6-31.6).As a result, the hitherto suspended thread with the handle hThread is resumed and executes the function ThreadProc carrying out the model calculations.In the case of a need to stop the simulation temporarily, the user may choose the button "Pause."Since the variable pause has the false value until this moment, the function SuspendThread is invoked, the pause reverts to true (lines 33.6-37.6),and all calculations in the function ThreadProc will stop.And, at last, if the user clicks a "Cancel" button, the dialog box procedure closes the dialog box by calling the EndDialog function (line 40.6) as it does in case of processing the WM CLOSE message (line 2.6).Next up is the case of processing another windows message WM HSCROLL that is sent if the user wants to change the simulation rate and moves the slider of the trackbar control to the left to increase the rate or to the right to slow it down.The simulation rate is controlled by the global variable sleep used further in the function ThreadProc.As mentioned above, sleep is assigned the value returned by the function SendMessage in lines 45.6-46.6.For information, this value is also shown at the edit box with an identifier IDC EDIT Sleep (lines 48.6-49.6).

Painting at the Modal Dialog Box: Processing the Message
WM PAINT.Drawing graphical objects at the modal dialog box is carried out when the system generates a message WM PAINT.It is done either implicitly when the window must be updated after sizing, moving, and so forth, or the application explicitly requests an update, invoking, for example, the function InvalidateRect.In our case, these key control actions are accomplished in the function ThreadProc (Algorithm 9, lines 36.9-38.9).
The typical processing of the WM PAINT starts with calling the function BeginPaint (Algorithm 6, line 51.6) which passes the painting information to the variable psPaint and returns the handle hDC to the device context.Then a number of functions are called to prepare the client area at the modal box to draw the interface elements and data (lines 52.6-55.6).When calling the function GetClientRect, a rectangle of the client area gets its coordinates; calling the function SetMap-Mode sets the mapping mode between page space and device space.In case of a parameter MM ISOTROPIC logical units are mapped to device (e.g., display) units with equally scaled axes.The function SetWindowExtEx sets the maximum values for the horizontal and vertical axes of the window in logical units, whereas the function SetViewportExtEx does the same in device units for the view port, that is, the device coordinate system.So, using both functions determines the scaling factor between the window and the view port.
Shown below (Algorithm 7, lines 2.7-5.7) is the fragment that draws the caption of the simulation application "A simulation model of the motion. ..." Next, in lines 6.7-50.7, the graphic functions are invoked such as Rectangle, Move-ToEx, and LineTo.These functions provide drawing a plot rectangle, axes, gridlines, and tick marks (see the screenshot in Figure 2(a)) necessary for displaying electrodes, field lines, particle's trajectories, and so on.The equipotential lines are drawn by calling the helper function DrawEquipotentialLines (lines 53.7-54.7).The principal painting, that is, drawing the particles' trajectories in real time, is carried out in the fragments placed between lines 2.8 and 14.8 in Algorithm 8.Both the global vector vpp and the Particle's data members, vectors history keeping the particles' states, are involved in the drawing of the trajectories () in the nested for loop.The processing of the message WM PAINT is completed with the fragment written in lines 17.8-21.8which displays the vpp's index of a particle whose data is to be output at the list box.The final statement of the processing is the call of EndPaint.
When the modal dialog box is being destroyed and is removed from the screen, the message WM DESTROY is sent to the function ModalDlgProc.The processing of this message includes the set of function calls of DeleteObject to delete the logical pens and the brush.

Description of the Modeless Dialog Box with the Console-Like
Output.The modeless dialog box can be shown even while the modal dialog box exists.Its purpose in this application is to display the numerical data describing the history of the state of a particle having been selected at the edit box with an identifier IDC EDIT ParticleModelessBox during the simulation (see Figure 2(b)).The output is carried out by the window function ModelessDlgProc in the console-like fashion and requires the vertical and horizontal scrolling to get more display space the modeless box can provide.That is, vertical and horizontal scroll bars are used in the modeless box.
The window function has the similar structure as in case of the modal dialog.It processes the messages in the switch selection statement that has eight branches, case labels for processing the messages WM INITDIALOG, WM VSCROLL, WM HSCROLL, WM PAINT, and so on.As the most of the procedure code contains the common, near standard fragments (see, e.g., Petzold [14], Chapter 4, section "Scroll Bars"), the listing of the procedure is placed to Appendix.However, a short description of the procedure is still given below.The resource script placed in the file script2.rcand associated with the modeless dialog contains styles WS VSCROLL and WS HSCROLL among other ones:

WS MINIMIZEBOX | WS OVERLAPPEDWINDOW
A number of variables and constants are defined in the window procedure.The constant NUMLINES is the maximal number of lines of text involved into the vertical scrolling.The variable tm of type TEXTMETRIC defines the font and rows parameters such as width and height of characters and space between characters.The static variables yPos and xPos are current positions of the scroll bar thumbs for vertical and horizontal scroll bars correspondingly, and yInc contains an incremental value for updating yPos.The messages WM VSCROLL and WM HSCROLL are sent to the window procedure if the scroll events occur either in the window's vertical scroll bar or in the horizontal scroll bar.Processing these messages made within nested switch statements (see Appendix, Algorithm 11, with the parameter notification code SB TOP, SB BOTTOM, and so forth results in updating yPos and xPos. When calling the function ScrollWindow in lines 68.11 and 96.11, the system invalidates the client area and generates the WM PAINT message.One more origin of a WM PAINT message is the thread procedure, where the function InvalidateRect is called (see Algorithm 9, line 38.9) after completing the nested for loop of one step integration of charged particles' motion.The code written in the WM PAINT selection brunch provides drawing the full set of text lines from the beginning of simulation to the current instant.Each text line outputs the time, position, velocity, and kinetic energy of a particle (i.e., the set of , , , V  , V  , V, ).

4.3.
Thread Procedure as the Physics Engine.The physics engine of the considered simulation model is placed in the function ThreadProc.The function includes the while loop (Algorithm 9, lines 4.9-40.9)providing the simulation until the last particle reaches the rightmost edge of the lens.The nested for loop (lines 5.9-25.9)carries out one step of integration for all particles which did not reach the right lens edge (i.e.,   < ).When a particle gets moved beyond this boundary, its data member stop is set to true; that is, the particle is excluded from simulation, and the variable counter keeping the number of excluded particles is incremented.The integration is carried out by the global function template RK4 4eq (lines 7.9-10.9).This function is a solver that processes one step of integration of the system of four ODEs by the Runge-Kutta method of the 4th order.The function template is defined in the header file rk4 eq.h and has the prototype as shown in Algorithm 1.
As we can see in the prototype shown in Algorithm 1, the function template imposes some restrictions on the content of the application class ; namely, the class needs to have five data members denoting both an independent variable (time ) and four state variables (, , V  , and V  ), a helper function (pointed as pPrecompute), and four member functions for computing the right-hand sides of the corresponding differential equations.The helper function makes the calculations common for right-hand sides, thus reducing the run time.The class Particle is designed so that these requirements are met.
After an integration step has been completed, the updated values of the state variables are pushed into the Particle's vector container history (lines 11.9-20.9).When the for loop is completed, the messages are generated to be passed both to the modal dialog box (including the list box) and to the modeless box.In particular, the functions InvalidateRect are invoked generating the messages WM PAINT.As a result, the numerical outputs are updated.When sending messages to the list box via calling the functions SendMessage in lines 31.9-34.9, the message LB ADDSTRING results in adding the new text line to the end of the list, and the message LB SETTOPINDEX ensures that the end of the list is visible.The call of the function Sleep allows setting the desired simulation rate kept in the variable sleep.Then the while loop starts the next iteration.was built in release configuration, so the program is fully optimized.The simulation is carried out on an AMD Athlon TM 64 × 2 Dual Core 2.40 GHz processor.With the aim of getting the CPU (central processing unit) time that elapses when executing the thread procedure ThreadProc (physics engine), the function clock defined in the header file <ctime> has been used.The time step of integration of the particle motion equations is ℎ  = 10 −9 sec.Three cases of simulation were considered, which have different sets of output facilities.The first case covers the program version with the graphical output only, the second one adds the text output to the modeless box, and the third one provides also the text output to the list box.In all cases the time delay is set to zero.The results obtained for the cases depending on the sizes of the particle ensemble are presented in Table 1.Each value is an average of a small sample formed as a result of 5-7 runs.The results show a high performance of the simulation model in the first (a) and the second (b) cases, even for large particle ensembles (10 4 particles).However, the third case (c), in which the text output is also displayed in the list box, is characterised by a much lower rate of computation and output.So this output mode should be applied only for applications not having the high requirements to the simulation speed.

Some Simulation Results and Discussion
1.9 DWORD CALLBACK ThreadProc(LPVOID lpdwThreadParam){ 2.9 const double ht = 1.0e-9; //time step of integration, s 3.9 static int count = 0; //number of particles with z > zMax 4.9 while(count < (signed)vpp.size()){5.9 for(int i = 0; i < (int)vpp.size();++i) 6.9 if Summing up this and other advantages of C++, such as greater flexibility through the wide use of pointers and templates and the availability of extensive libraries (Win32 API, STL, Boost, etc.), we have every reason to consider C++ as a premier object-oriented language for scientific applications.

Conclusion
The case study of the simulation model considered in the paper provides some general guidelines on how to develop a desktop framework that would have the minimal set of means being sufficient for organizing the simulation process as such,

Figure 1 :
Figure 1: The block diagram showing relationships between classes, objects, functions, and global variables being included in the simulation model.

3. 1 .Figure 2 :
Figure 2: Screenshot of the graphic user interface (a) and a magnified view of the part of modal dialog box with control items (b).
Modal Dialog Procedure.Next we describe the modal dialog box window procedure ModalDl-gProc; the definition of that is listed in Algorithms 5-8.The function has the following heading: int CALLBACK ModalDlgProc(HWND hwnd, UINT message, WPARAM wParam, LPARAM lParam ).
(v) An analogous procedure accomplished concerning the trackbar control which provides changing the rate of simulation: the functions SendMessage in the lines 48.5-51.5 define the range and set the slider to the initial position between two tick marks; another SendMessage function gets the slider's logical position that is assigned to the global variable sleep (lines 52.5-53.5).(vi)Statements in lines 42.5-47.5 and in lines 54.5-55.5 setting the initial values for all edit windows.
4.1.1.Processing the Control Related Messages WM COMMAND.Another important case we need to focus on is processing the message WM COMMAND which is sent when the user clicks a button at the modal dialog box.If the user chooses a button "Create a new particle," the function GetDlgItemText retrieves strings from the corresponding edit windows with control identifiers IDC EDIT W0 and IDC EDIT kVy to the text buffer (Algorithm 6, lines 8.6 and 10.6).Then the strings are converted into numerical values of the variables W0 and kVy (lines 9.6 and 11.6) and are passed to the Particle's constructor.The latter creates a Particle's object that will be kept in the global STL vector container vpp: 12.6 vpp.push back(Particle(W0, kVy)); The considered object-oriented simulation model of particles motion in electric field combined with Win32 dialog-based graphical interface and framework has been tested within Visual C++ 2010 Express environment for different cases to estimate the performance of the simulation.The project template<class T> inline void RK4 4eq(T* p, //pointer to an object of the class T double h, //time step of integration (ht in the class Particle) double T::*pX0, // pointer to the independent variable -time  double T::*pX1, // pointer to the state variable -velocity V  Algorithm 2: The header file particle.hdefining the class Particle.

Table 1 :
The elapsed CPU time required for simulating the motion of an ensemble of particles in the thread function ThreadProc (no sleep time).