88 * @details
99 * This example demonstrates how to create custom field initialization patterns
1010 * using simple structs and evaluation functions. We show three physical patterns:
11- *
11+ *
1212 * 1. **Lamb-Oseen Vortex** - Rotating fluid vortex with viscous core
1313 * 2. **Gaussian Bump** - Localized concentration or temperature field
1414 * 3. **Checkerboard** - Periodic alternating pattern
2727 * This simplified example focuses on the pattern concept itself.
2828 */
2929
30- #include < openpfc/openpfc.hpp>
3130#include < cmath>
3231#include < iostream>
32+ #include < openpfc/openpfc.hpp>
3333
3434using namespace pfc ;
3535
@@ -49,48 +49,48 @@ namespace my_project {
4949
5050/* *
5151 * @brief Lamb-Oseen vortex pattern
52- *
52+ *
5353 * Models a rotating vortex with viscous core, common in fluid dynamics.
5454 * The tangential velocity follows: v_θ(r) = (Γ/2πr)[1 - exp(-r²/r_c²)]
55- *
55+ *
5656 * @see https://en.wikipedia.org/wiki/Lamb%E2%80%93Oseen_vortex
5757 */
5858struct VortexPattern {
59- Real3 m_center; // /< Vortex center (x, y, z)
60- double m_strength; // /< Circulation Γ
61- double m_core_radius; // /< Core radius r_c
62-
63- VortexPattern (Real3 center, double strength, double core_radius)
64- : m_center(center), m_strength(strength), m_core_radius(core_radius) {}
59+ Real3 m_center; // /< Vortex center (x, y, z)
60+ double m_strength; // /< Circulation Γ
61+ double m_core_radius; // /< Core radius r_c
62+
63+ VortexPattern (Real3 center, double strength, double core_radius)
64+ : m_center(center), m_strength(strength), m_core_radius(core_radius) {}
6565};
6666
6767/* *
6868 * @brief 3D Gaussian bump
69- *
69+ *
7070 * Models localized concentration, temperature, or density field.
7171 * φ(r) = A * exp(-r²/(2σ²))
7272 */
7373struct GaussianBump {
74- Real3 m_center; // /< Peak location
75- double m_amplitude; // /< Peak height A
76- double m_width; // /< Standard deviation σ
77-
78- GaussianBump (Real3 center, double amplitude, double width)
79- : m_center(center), m_amplitude(amplitude), m_width(width) {}
74+ Real3 m_center; // /< Peak location
75+ double m_amplitude; // /< Peak height A
76+ double m_width; // /< Standard deviation σ
77+
78+ GaussianBump (Real3 center, double amplitude, double width)
79+ : m_center(center), m_amplitude(amplitude), m_width(width) {}
8080};
8181
8282/* *
8383 * @brief 3D checkerboard pattern
84- *
84+ *
8585 * Periodic alternating values, useful for testing and validation.
8686 */
8787struct CheckerboardPattern {
88- double m_value_high; // /< Value in "white" cells
89- double m_value_low; // /< Value in "black" cells
90- Real3 m_period; // /< Period in each direction
91-
92- CheckerboardPattern (double high, double low, Real3 period)
93- : m_value_high(high), m_value_low(low), m_period(period) {}
88+ double m_value_high; // /< Value in "white" cells
89+ double m_value_low; // /< Value in "black" cells
90+ Real3 m_period; // /< Period in each direction
91+
92+ CheckerboardPattern (double high, double low, Real3 period)
93+ : m_value_high(high), m_value_low(low), m_period(period) {}
9494};
9595
9696} // namespace my_project
@@ -105,22 +105,22 @@ struct CheckerboardPattern {
105105 * @param pos Physical position to evaluate at
106106 * @return Tangential velocity at position
107107 */
108- double evaluate_vortex (const my_project::VortexPattern& pattern, const Real3& pos) {
109- // Distance from vortex center in x-y plane
110- double dx = pos[0 ] - pattern.m_center [0 ];
111- double dy = pos[1 ] - pattern.m_center [1 ];
112- double r = std::sqrt (dx * dx + dy * dy);
113-
114- // Lamb-Oseen vortex profile
115- double r_c_sq = pattern.m_core_radius * pattern.m_core_radius ;
116- double value = 0.0 ;
117-
118- if (r > 1e-10 ) { // Avoid division by zero at center
119- value = (pattern. m_strength / ( 2.0 * M_PI * r)) *
120- (1.0 - std::exp (-r * r / r_c_sq));
121- }
122-
123- return value;
108+ double evaluate_vortex (const my_project::VortexPattern & pattern, const Real3 & pos) {
109+ // Distance from vortex center in x-y plane
110+ double dx = pos[0 ] - pattern.m_center [0 ];
111+ double dy = pos[1 ] - pattern.m_center [1 ];
112+ double r = std::sqrt (dx * dx + dy * dy);
113+
114+ // Lamb-Oseen vortex profile
115+ double r_c_sq = pattern.m_core_radius * pattern.m_core_radius ;
116+ double value = 0.0 ;
117+
118+ if (r > 1e-10 ) { // Avoid division by zero at center
119+ value =
120+ (pattern. m_strength / ( 2.0 * M_PI * r)) * (1.0 - std::exp (-r * r / r_c_sq));
121+ }
122+
123+ return value;
124124}
125125
126126/* *
@@ -129,16 +129,16 @@ double evaluate_vortex(const my_project::VortexPattern& pattern, const Real3& po
129129 * @param pos Physical position to evaluate at
130130 * @return Field value at position
131131 */
132- double evaluate_gaussian (const my_project::GaussianBump& pattern, const Real3& pos) {
133- // Distance from center
134- double dx = pos[0 ] - pattern.m_center [0 ];
135- double dy = pos[1 ] - pattern.m_center [1 ];
136- double dz = pos[2 ] - pattern.m_center [2 ];
137- double dist_sq = dx* dx + dy* dy + dz* dz;
138-
139- // Gaussian: φ = A * exp(-dist² / (2σ²))
140- double sigma_sq = pattern.m_width * pattern.m_width ;
141- return pattern.m_amplitude * std::exp (-dist_sq / (2.0 * sigma_sq));
132+ double evaluate_gaussian (const my_project::GaussianBump & pattern, const Real3 & pos) {
133+ // Distance from center
134+ double dx = pos[0 ] - pattern.m_center [0 ];
135+ double dy = pos[1 ] - pattern.m_center [1 ];
136+ double dz = pos[2 ] - pattern.m_center [2 ];
137+ double dist_sq = dx * dx + dy * dy + dz * dz;
138+
139+ // Gaussian: φ = A * exp(-dist² / (2σ²))
140+ double sigma_sq = pattern.m_width * pattern.m_width ;
141+ return pattern.m_amplitude * std::exp (-dist_sq / (2.0 * sigma_sq));
142142}
143143
144144/* *
@@ -147,90 +147,91 @@ double evaluate_gaussian(const my_project::GaussianBump& pattern, const Real3& p
147147 * @param pos Physical position to evaluate at
148148 * @return High or low value depending on position
149149 */
150- double evaluate_checkerboard (const my_project::CheckerboardPattern& pattern, const Real3& pos) {
151- // Determine which cell of the checkerboard
152- int cell_i = static_cast <int >(std::floor (pos[0 ] / pattern.m_period [0 ]));
153- int cell_j = static_cast <int >(std::floor (pos[1 ] / pattern.m_period [1 ]));
154- int cell_k = static_cast <int >(std::floor (pos[2 ] / pattern.m_period [2 ]));
155-
156- // Checkerboard: alternate based on sum of cell indices
157- int sum = cell_i + cell_j + cell_k;
158- return (sum % 2 == 0 ) ? pattern.m_value_high : pattern.m_value_low ;
150+ double evaluate_checkerboard (const my_project::CheckerboardPattern &pattern,
151+ const Real3 &pos) {
152+ // Determine which cell of the checkerboard
153+ int cell_i = static_cast <int >(std::floor (pos[0 ] / pattern.m_period [0 ]));
154+ int cell_j = static_cast <int >(std::floor (pos[1 ] / pattern.m_period [1 ]));
155+ int cell_k = static_cast <int >(std::floor (pos[2 ] / pattern.m_period [2 ]));
156+
157+ // Checkerboard: alternate based on sum of cell indices
158+ int sum = cell_i + cell_j + cell_k;
159+ return (sum % 2 == 0 ) ? pattern.m_value_high : pattern.m_value_low ;
159160}
160161
161162// ============================================================================
162163// Part 3: Usage Examples
163164// ============================================================================
164165
165166void example_vortex_pattern () {
166- std::cout << " === Example 1: Vortex Pattern ===\n\n " ;
167-
168- // Define vortex at origin with strength 5.0 and core radius 2.0
169- my_project::VortexPattern vortex ({0.0 , 0.0 , 0.0 }, 5.0 , 2.0 );
170-
171- // Evaluate at several radial distances
172- std::cout << " Vortex tangential velocity profile:\n " ;
173- for (double r = 0.0 ; r <= 10.0 ; r += 2.0 ) {
174- Real3 pos{r, 0.0 , 0.0 }; // Along x-axis
175- double velocity = evaluate_vortex (vortex, pos);
176- std::cout << " r = " << r << " : v_θ = " << velocity << " \n " ;
177- }
178- std::cout << " \n " ;
167+ std::cout << " === Example 1: Vortex Pattern ===\n\n " ;
168+
169+ // Define vortex at origin with strength 5.0 and core radius 2.0
170+ my_project::VortexPattern vortex ({0.0 , 0.0 , 0.0 }, 5.0 , 2.0 );
171+
172+ // Evaluate at several radial distances
173+ std::cout << " Vortex tangential velocity profile:\n " ;
174+ for (double r = 0.0 ; r <= 10.0 ; r += 2.0 ) {
175+ Real3 pos{r, 0.0 , 0.0 }; // Along x-axis
176+ double velocity = evaluate_vortex (vortex, pos);
177+ std::cout << " r = " << r << " : v_θ = " << velocity << " \n " ;
178+ }
179+ std::cout << " \n " ;
179180}
180181
181182void example_gaussian_bump () {
182- std::cout << " === Example 2: Gaussian Bump ===\n\n " ;
183-
184- // Create Gaussian peak at origin with amplitude 1.0, width 1.5
185- my_project::GaussianBump bump ({0.0 , 0.0 , 0.0 }, 1.0 , 1.5 );
186-
187- // Evaluate along a line
188- std::cout << " Gaussian profile:\n " ;
189- for (double x = 0.0 ; x <= 5.0 ; x += 1.0 ) {
190- Real3 pos{x, 0.0 , 0.0 };
191- double value = evaluate_gaussian (bump, pos);
192- std::cout << " x = " << x << " : φ = " << value << " \n " ;
193- }
194- std::cout << " \n " ;
183+ std::cout << " === Example 2: Gaussian Bump ===\n\n " ;
184+
185+ // Create Gaussian peak at origin with amplitude 1.0, width 1.5
186+ my_project::GaussianBump bump ({0.0 , 0.0 , 0.0 }, 1.0 , 1.5 );
187+
188+ // Evaluate along a line
189+ std::cout << " Gaussian profile:\n " ;
190+ for (double x = 0.0 ; x <= 5.0 ; x += 1.0 ) {
191+ Real3 pos{x, 0.0 , 0.0 };
192+ double value = evaluate_gaussian (bump, pos);
193+ std::cout << " x = " << x << " : φ = " << value << " \n " ;
194+ }
195+ std::cout << " \n " ;
195196}
196197
197198void example_checkerboard () {
198- std::cout << " === Example 3: Checkerboard Pattern ===\n\n " ;
199-
200- // Create checkerboard with period 2.0 in each direction
201- my_project::CheckerboardPattern checker (1.0 , -1.0 , {2.0 , 2.0 , 2.0 });
202-
203- // Sample a 2D slice (z=0 plane)
204- std::cout << " Checkerboard pattern (z=0 plane):\n " ;
205- for (int j = 0 ; j < 4 ; ++j) {
206- std::cout << " " ;
207- for (int i = 0 ; i < 4 ; ++i) {
208- Real3 pos{i * 1.0 , j * 1.0 , 0.0 };
209- double value = evaluate_checkerboard (checker, pos);
210- std::cout << (value > 0 ? " +" : " -" ) << " " ;
211- }
212- std::cout << " \n " ;
199+ std::cout << " === Example 3: Checkerboard Pattern ===\n\n " ;
200+
201+ // Create checkerboard with period 2.0 in each direction
202+ my_project::CheckerboardPattern checker (1.0 , -1.0 , {2.0 , 2.0 , 2.0 });
203+
204+ // Sample a 2D slice (z=0 plane)
205+ std::cout << " Checkerboard pattern (z=0 plane):\n " ;
206+ for (int j = 0 ; j < 4 ; ++j) {
207+ std::cout << " " ;
208+ for (int i = 0 ; i < 4 ; ++i) {
209+ Real3 pos{i * 1.0 , j * 1.0 , 0.0 };
210+ double value = evaluate_checkerboard (checker, pos);
211+ std::cout << (value > 0 ? " +" : " -" ) << " " ;
213212 }
214213 std::cout << " \n " ;
214+ }
215+ std::cout << " \n " ;
215216}
216217
217218int main () {
218- std::cout << " \n " ;
219- std::cout << " ╔════════════════════════════════════════════════════════╗\n " ;
220- std::cout << " ║ OpenPFC: Custom Field Initialization Patterns Example ║\n " ;
221- std::cout << " ╚════════════════════════════════════════════════════════╝\n " ;
222- std::cout << " \n " ;
223-
224- example_vortex_pattern ();
225- example_gaussian_bump ();
226- example_checkerboard ();
227-
228- std::cout << " ✅ All examples completed successfully!\n " ;
229- std::cout << " \n " ;
230- std::cout << " 📖 For full integration with DiscreteField/World using ADL, see:\n " ;
231- std::cout << " - docs/extending_openpfc/adl_extension_patterns.md\n " ;
232- std::cout << " - examples/17_custom_coordinate_system.cpp\n " ;
233- std::cout << " \n " ;
234-
235- return 0 ;
219+ std::cout << " \n " ;
220+ std::cout << " ╔════════════════════════════════════════════════════════╗\n " ;
221+ std::cout << " ║ OpenPFC: Custom Field Initialization Patterns Example ║\n " ;
222+ std::cout << " ╚════════════════════════════════════════════════════════╝\n " ;
223+ std::cout << " \n " ;
224+
225+ example_vortex_pattern ();
226+ example_gaussian_bump ();
227+ example_checkerboard ();
228+
229+ std::cout << " ✅ All examples completed successfully!\n " ;
230+ std::cout << " \n " ;
231+ std::cout << " 📖 For full integration with DiscreteField/World using ADL, see:\n " ;
232+ std::cout << " - docs/extending_openpfc/adl_extension_patterns.md\n " ;
233+ std::cout << " - examples/17_custom_coordinate_system.cpp\n " ;
234+ std::cout << " \n " ;
235+
236+ return 0 ;
236237}
0 commit comments