1 // Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at 2 // the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights 3 // reserved. See files LICENSE and NOTICE for details. 4 // 5 // This file is part of CEED, a collection of benchmarks, miniapps, software 6 // libraries and APIs for efficient high-order finite element and spectral 7 // element discretizations for exascale applications. For more information and 8 // source code availability see http://github.com/ceed. 9 // 10 // The CEED research is supported by the Exascale Computing Project 17-SC-20-SC, 11 // a collaborative effort of two U.S. Department of Energy organizations (Office 12 // of Science and the National Nuclear Security Administration) responsible for 13 // the planning and preparation of a capable exascale ecosystem, including 14 // software, applications, hardware, advanced system engineering and early 15 // testbed platforms, in support of the nation's exascale computing imperative 16 17 //! A Ceed QFunction represents the spatial terms of the point-wise functions 18 //! describing the physics at the quadrature points. 19 20 use std::pin::Pin; 21 22 use crate::prelude::*; 23 24 pub type QFunctionInputs<'a> = [&'a [crate::Scalar]; MAX_QFUNCTION_FIELDS]; 25 pub type QFunctionOutputs<'a> = [&'a mut [crate::Scalar]; MAX_QFUNCTION_FIELDS]; 26 27 // ----------------------------------------------------------------------------- 28 // CeedQFunction option 29 // ----------------------------------------------------------------------------- 30 pub enum QFunctionOpt<'a> { 31 SomeQFunction(&'a QFunction<'a>), 32 SomeQFunctionByName(&'a QFunctionByName<'a>), 33 None, 34 } 35 36 /// Construct a QFunctionOpt reference from a QFunction reference 37 impl<'a> From<&'a QFunction<'_>> for QFunctionOpt<'a> { 38 fn from(qfunc: &'a QFunction) -> Self { 39 debug_assert!(qfunc.qf_core.ptr != unsafe { bind_ceed::CEED_QFUNCTION_NONE }); 40 Self::SomeQFunction(qfunc) 41 } 42 } 43 44 /// Construct a QFunctionOpt reference from a QFunction by Name reference 45 impl<'a> From<&'a QFunctionByName<'_>> for QFunctionOpt<'a> { 46 fn from(qfunc: &'a QFunctionByName) -> Self { 47 debug_assert!(qfunc.qf_core.ptr != unsafe { bind_ceed::CEED_QFUNCTION_NONE }); 48 Self::SomeQFunctionByName(qfunc) 49 } 50 } 51 52 impl<'a> QFunctionOpt<'a> { 53 /// Transform a Rust libCEED QFunctionOpt into C libCEED CeedQFunction 54 pub(crate) fn to_raw(self) -> bind_ceed::CeedQFunction { 55 match self { 56 Self::SomeQFunction(qfunc) => qfunc.qf_core.ptr, 57 Self::SomeQFunctionByName(qfunc) => qfunc.qf_core.ptr, 58 Self::None => unsafe { bind_ceed::CEED_QFUNCTION_NONE }, 59 } 60 } 61 62 /// Check if a QFunctionOpt is Some 63 /// 64 /// ``` 65 /// # use libceed::prelude::*; 66 /// # fn main() -> Result<(), libceed::CeedError> { 67 /// # let ceed = libceed::Ceed::default_init(); 68 /// let mut user_f = |[u, weights, ..]: QFunctionInputs, [v, ..]: QFunctionOutputs| { 69 /// // Iterate over quadrature points 70 /// v.iter_mut() 71 /// .zip(u.iter().zip(weights.iter())) 72 /// .for_each(|(v, (u, w))| *v = u * w); 73 /// 74 /// // Return clean error code 75 /// 0 76 /// }; 77 /// 78 /// let qf = ceed 79 /// .q_function_interior(1, Box::new(user_f))? 80 /// .input("u", 1, EvalMode::Interp)? 81 /// .input("weights", 1, EvalMode::Weight)? 82 /// .output("v", 1, EvalMode::Interp)?; 83 /// let qf_opt = QFunctionOpt::from(&qf); 84 /// assert!(qf_opt.is_some(), "Incorrect QFunctionOpt"); 85 /// 86 /// let qf = ceed.q_function_interior_by_name("Mass1DBuild")?; 87 /// let qf_opt = QFunctionOpt::from(&qf); 88 /// assert!(qf_opt.is_some(), "Incorrect QFunctionOpt"); 89 /// 90 /// let qf_opt = QFunctionOpt::None; 91 /// assert!(!qf_opt.is_some(), "Incorrect QFunctionOpt"); 92 /// # Ok(()) 93 /// # } 94 /// ``` 95 pub fn is_some(&self) -> bool { 96 match self { 97 Self::SomeQFunction(_) => true, 98 Self::SomeQFunctionByName(_) => true, 99 Self::None => false, 100 } 101 } 102 103 /// Check if a QFunctionOpt is SomeQFunction 104 /// 105 /// ``` 106 /// # use libceed::prelude::*; 107 /// # fn main() -> Result<(), libceed::CeedError> { 108 /// # let ceed = libceed::Ceed::default_init(); 109 /// let mut user_f = |[u, weights, ..]: QFunctionInputs, [v, ..]: QFunctionOutputs| { 110 /// // Iterate over quadrature points 111 /// v.iter_mut() 112 /// .zip(u.iter().zip(weights.iter())) 113 /// .for_each(|(v, (u, w))| *v = u * w); 114 /// 115 /// // Return clean error code 116 /// 0 117 /// }; 118 /// 119 /// let qf = ceed 120 /// .q_function_interior(1, Box::new(user_f))? 121 /// .input("u", 1, EvalMode::Interp)? 122 /// .input("weights", 1, EvalMode::Weight)? 123 /// .output("v", 1, EvalMode::Interp)?; 124 /// let qf_opt = QFunctionOpt::from(&qf); 125 /// assert!(qf_opt.is_some_q_function(), "Incorrect QFunctionOpt"); 126 /// 127 /// let qf = ceed.q_function_interior_by_name("Mass1DBuild")?; 128 /// let qf_opt = QFunctionOpt::from(&qf); 129 /// assert!(!qf_opt.is_some_q_function(), "Incorrect QFunctionOpt"); 130 /// 131 /// let qf_opt = QFunctionOpt::None; 132 /// assert!(!qf_opt.is_some_q_function(), "Incorrect QFunctionOpt"); 133 /// # Ok(()) 134 /// # } 135 /// ``` 136 pub fn is_some_q_function(&self) -> bool { 137 match self { 138 Self::SomeQFunction(_) => true, 139 Self::SomeQFunctionByName(_) => false, 140 Self::None => false, 141 } 142 } 143 144 /// Check if a QFunctionOpt is SomeQFunctionByName 145 /// 146 /// ``` 147 /// # use libceed::prelude::*; 148 /// # fn main() -> Result<(), libceed::CeedError> { 149 /// # let ceed = libceed::Ceed::default_init(); 150 /// let mut user_f = |[u, weights, ..]: QFunctionInputs, [v, ..]: QFunctionOutputs| { 151 /// // Iterate over quadrature points 152 /// v.iter_mut() 153 /// .zip(u.iter().zip(weights.iter())) 154 /// .for_each(|(v, (u, w))| *v = u * w); 155 /// 156 /// // Return clean error code 157 /// 0 158 /// }; 159 /// 160 /// let qf = ceed 161 /// .q_function_interior(1, Box::new(user_f))? 162 /// .input("u", 1, EvalMode::Interp)? 163 /// .input("weights", 1, EvalMode::Weight)? 164 /// .output("v", 1, EvalMode::Interp)?; 165 /// let qf_opt = QFunctionOpt::from(&qf); 166 /// assert!( 167 /// !qf_opt.is_some_q_function_by_name(), 168 /// "Incorrect QFunctionOpt" 169 /// ); 170 /// 171 /// let qf = ceed.q_function_interior_by_name("Mass1DBuild")?; 172 /// let qf_opt = QFunctionOpt::from(&qf); 173 /// assert!( 174 /// qf_opt.is_some_q_function_by_name(), 175 /// "Incorrect QFunctionOpt" 176 /// ); 177 /// 178 /// let qf_opt = QFunctionOpt::None; 179 /// assert!( 180 /// !qf_opt.is_some_q_function_by_name(), 181 /// "Incorrect QFunctionOpt" 182 /// ); 183 /// # Ok(()) 184 /// # } 185 /// ``` 186 pub fn is_some_q_function_by_name(&self) -> bool { 187 match self { 188 Self::SomeQFunction(_) => false, 189 Self::SomeQFunctionByName(_) => true, 190 Self::None => false, 191 } 192 } 193 194 /// Check if a QFunctionOpt is None 195 /// 196 /// ``` 197 /// # use libceed::prelude::*; 198 /// # fn main() -> Result<(), libceed::CeedError> { 199 /// # let ceed = libceed::Ceed::default_init(); 200 /// let mut user_f = |[u, weights, ..]: QFunctionInputs, [v, ..]: QFunctionOutputs| { 201 /// // Iterate over quadrature points 202 /// v.iter_mut() 203 /// .zip(u.iter().zip(weights.iter())) 204 /// .for_each(|(v, (u, w))| *v = u * w); 205 /// 206 /// // Return clean error code 207 /// 0 208 /// }; 209 /// 210 /// let qf = ceed 211 /// .q_function_interior(1, Box::new(user_f))? 212 /// .input("u", 1, EvalMode::Interp)? 213 /// .input("weights", 1, EvalMode::Weight)? 214 /// .output("v", 1, EvalMode::Interp)?; 215 /// let qf_opt = QFunctionOpt::from(&qf); 216 /// assert!(!qf_opt.is_none(), "Incorrect QFunctionOpt"); 217 /// 218 /// let qf = ceed.q_function_interior_by_name("Mass1DBuild")?; 219 /// let qf_opt = QFunctionOpt::from(&qf); 220 /// assert!(!qf_opt.is_none(), "Incorrect QFunctionOpt"); 221 /// 222 /// let qf_opt = QFunctionOpt::None; 223 /// assert!(qf_opt.is_none(), "Incorrect QFunctionOpt"); 224 /// # Ok(()) 225 /// # } 226 /// ``` 227 pub fn is_none(&self) -> bool { 228 match self { 229 Self::SomeQFunction(_) => false, 230 Self::SomeQFunctionByName(_) => false, 231 Self::None => true, 232 } 233 } 234 } 235 236 // ----------------------------------------------------------------------------- 237 // CeedQFunction context wrapper 238 // ----------------------------------------------------------------------------- 239 #[derive(Debug)] 240 pub(crate) struct QFunctionCore<'a> { 241 ceed: &'a crate::Ceed, 242 ptr: bind_ceed::CeedQFunction, 243 } 244 245 struct QFunctionTrampolineData { 246 number_inputs: usize, 247 number_outputs: usize, 248 input_sizes: [usize; MAX_QFUNCTION_FIELDS], 249 output_sizes: [usize; MAX_QFUNCTION_FIELDS], 250 user_f: Box<QFunctionUserClosure>, 251 } 252 253 pub struct QFunction<'a> { 254 qf_core: QFunctionCore<'a>, 255 qf_ctx_ptr: bind_ceed::CeedQFunctionContext, 256 trampoline_data: Pin<Box<QFunctionTrampolineData>>, 257 } 258 259 #[derive(Debug)] 260 pub struct QFunctionByName<'a> { 261 qf_core: QFunctionCore<'a>, 262 } 263 264 // ----------------------------------------------------------------------------- 265 // Destructor 266 // ----------------------------------------------------------------------------- 267 impl<'a> Drop for QFunctionCore<'a> { 268 fn drop(&mut self) { 269 unsafe { 270 if self.ptr != bind_ceed::CEED_QFUNCTION_NONE { 271 bind_ceed::CeedQFunctionDestroy(&mut self.ptr); 272 } 273 } 274 } 275 } 276 277 impl<'a> Drop for QFunction<'a> { 278 fn drop(&mut self) { 279 unsafe { 280 bind_ceed::CeedQFunctionContextDestroy(&mut self.qf_ctx_ptr); 281 } 282 } 283 } 284 285 // ----------------------------------------------------------------------------- 286 // Display 287 // ----------------------------------------------------------------------------- 288 impl<'a> fmt::Display for QFunctionCore<'a> { 289 fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { 290 let mut ptr = std::ptr::null_mut(); 291 let mut sizeloc = crate::MAX_BUFFER_LENGTH; 292 let cstring = unsafe { 293 let file = bind_ceed::open_memstream(&mut ptr, &mut sizeloc); 294 bind_ceed::CeedQFunctionView(self.ptr, file); 295 bind_ceed::fclose(file); 296 CString::from_raw(ptr) 297 }; 298 cstring.to_string_lossy().fmt(f) 299 } 300 } 301 /// View a QFunction 302 /// 303 /// ``` 304 /// # use libceed::prelude::*; 305 /// # fn main() -> Result<(), libceed::CeedError> { 306 /// # let ceed = libceed::Ceed::default_init(); 307 /// let mut user_f = |[u, weights, ..]: QFunctionInputs, [v, ..]: QFunctionOutputs| { 308 /// // Iterate over quadrature points 309 /// v.iter_mut() 310 /// .zip(u.iter().zip(weights.iter())) 311 /// .for_each(|(v, (u, w))| *v = u * w); 312 /// 313 /// // Return clean error code 314 /// 0 315 /// }; 316 /// 317 /// let qf = ceed 318 /// .q_function_interior(1, Box::new(user_f))? 319 /// .input("u", 1, EvalMode::Interp)? 320 /// .input("weights", 1, EvalMode::Weight)? 321 /// .output("v", 1, EvalMode::Interp)?; 322 /// 323 /// println!("{}", qf); 324 /// # Ok(()) 325 /// # } 326 /// ``` 327 impl<'a> fmt::Display for QFunction<'a> { 328 fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { 329 self.qf_core.fmt(f) 330 } 331 } 332 333 /// View a QFunction by Name 334 /// 335 /// ``` 336 /// # use libceed::prelude::*; 337 /// # fn main() -> Result<(), libceed::CeedError> { 338 /// # let ceed = libceed::Ceed::default_init(); 339 /// let qf = ceed.q_function_interior_by_name("Mass1DBuild")?; 340 /// println!("{}", qf); 341 /// # Ok(()) 342 /// # } 343 /// ``` 344 impl<'a> fmt::Display for QFunctionByName<'a> { 345 fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { 346 self.qf_core.fmt(f) 347 } 348 } 349 350 // ----------------------------------------------------------------------------- 351 // Core functionality 352 // ----------------------------------------------------------------------------- 353 impl<'a> QFunctionCore<'a> { 354 // Common implementation 355 pub fn apply(&self, Q: usize, u: &[Vector], v: &[Vector]) -> crate::Result<i32> { 356 let mut u_c = [std::ptr::null_mut(); MAX_QFUNCTION_FIELDS]; 357 for i in 0..std::cmp::min(MAX_QFUNCTION_FIELDS, u.len()) { 358 u_c[i] = u[i].ptr; 359 } 360 let mut v_c = [std::ptr::null_mut(); MAX_QFUNCTION_FIELDS]; 361 for i in 0..std::cmp::min(MAX_QFUNCTION_FIELDS, v.len()) { 362 v_c[i] = v[i].ptr; 363 } 364 let Q = i32::try_from(Q).unwrap(); 365 let ierr = unsafe { 366 bind_ceed::CeedQFunctionApply(self.ptr, Q, u_c.as_mut_ptr(), v_c.as_mut_ptr()) 367 }; 368 self.ceed.check_error(ierr) 369 } 370 } 371 372 // ----------------------------------------------------------------------------- 373 // User QFunction Closure 374 // ----------------------------------------------------------------------------- 375 pub type QFunctionUserClosure = dyn FnMut( 376 [&[crate::Scalar]; MAX_QFUNCTION_FIELDS], 377 [&mut [crate::Scalar]; MAX_QFUNCTION_FIELDS], 378 ) -> i32; 379 380 macro_rules! mut_max_fields { 381 ($e:expr) => { 382 [ 383 $e, $e, $e, $e, $e, $e, $e, $e, $e, $e, $e, $e, $e, $e, $e, $e, 384 ] 385 }; 386 } 387 unsafe extern "C" fn trampoline( 388 ctx: *mut ::std::os::raw::c_void, 389 q: bind_ceed::CeedInt, 390 inputs: *const *const bind_ceed::CeedScalar, 391 outputs: *const *mut bind_ceed::CeedScalar, 392 ) -> ::std::os::raw::c_int { 393 let trampoline_data: Pin<&mut QFunctionTrampolineData> = std::mem::transmute(ctx); 394 395 // Inputs 396 let inputs_slice: &[*const bind_ceed::CeedScalar] = 397 std::slice::from_raw_parts(inputs, MAX_QFUNCTION_FIELDS); 398 let mut inputs_array: [&[crate::Scalar]; MAX_QFUNCTION_FIELDS] = [&[0.0]; MAX_QFUNCTION_FIELDS]; 399 inputs_slice 400 .iter() 401 .enumerate() 402 .map(|(i, &x)| { 403 std::slice::from_raw_parts(x, trampoline_data.input_sizes[i] * q as usize) 404 as &[crate::Scalar] 405 }) 406 .zip(inputs_array.iter_mut()) 407 .for_each(|(x, a)| *a = x); 408 409 // Outputs 410 let outputs_slice: &[*mut bind_ceed::CeedScalar] = 411 std::slice::from_raw_parts(outputs, MAX_QFUNCTION_FIELDS); 412 let mut outputs_array: [&mut [crate::Scalar]; MAX_QFUNCTION_FIELDS] = 413 mut_max_fields!(&mut [0.0]); 414 outputs_slice 415 .iter() 416 .enumerate() 417 .map(|(i, &x)| { 418 std::slice::from_raw_parts_mut(x, trampoline_data.output_sizes[i] * q as usize) 419 as &mut [crate::Scalar] 420 }) 421 .zip(outputs_array.iter_mut()) 422 .for_each(|(x, a)| *a = x); 423 424 // User closure 425 (trampoline_data.get_unchecked_mut().user_f)(inputs_array, outputs_array) 426 } 427 428 // ----------------------------------------------------------------------------- 429 // QFunction 430 // ----------------------------------------------------------------------------- 431 impl<'a> QFunction<'a> { 432 // Constructor 433 pub fn create( 434 ceed: &'a crate::Ceed, 435 vlength: usize, 436 user_f: Box<QFunctionUserClosure>, 437 ) -> crate::Result<Self> { 438 let source_c = CString::new("").expect("CString::new failed"); 439 let mut ptr = std::ptr::null_mut(); 440 441 // Context for closure 442 let number_inputs = 0; 443 let number_outputs = 0; 444 let input_sizes = [0; MAX_QFUNCTION_FIELDS]; 445 let output_sizes = [0; MAX_QFUNCTION_FIELDS]; 446 let trampoline_data = unsafe { 447 Pin::new_unchecked(Box::new(QFunctionTrampolineData { 448 number_inputs, 449 number_outputs, 450 input_sizes, 451 output_sizes, 452 user_f, 453 })) 454 }; 455 456 // Create QFunction 457 let vlength = i32::try_from(vlength).unwrap(); 458 let mut ierr = unsafe { 459 bind_ceed::CeedQFunctionCreateInterior( 460 ceed.ptr, 461 vlength, 462 Some(trampoline), 463 source_c.as_ptr(), 464 &mut ptr, 465 ) 466 }; 467 ceed.check_error(ierr)?; 468 469 // Set closure 470 let mut qf_ctx_ptr = std::ptr::null_mut(); 471 ierr = unsafe { bind_ceed::CeedQFunctionContextCreate(ceed.ptr, &mut qf_ctx_ptr) }; 472 ceed.check_error(ierr)?; 473 ierr = unsafe { 474 bind_ceed::CeedQFunctionContextSetData( 475 qf_ctx_ptr, 476 crate::MemType::Host as bind_ceed::CeedMemType, 477 crate::CopyMode::UsePointer as bind_ceed::CeedCopyMode, 478 std::mem::size_of::<QFunctionTrampolineData>() as u64, 479 std::mem::transmute(trampoline_data.as_ref()), 480 ) 481 }; 482 ceed.check_error(ierr)?; 483 ierr = unsafe { bind_ceed::CeedQFunctionSetContext(ptr, qf_ctx_ptr) }; 484 ceed.check_error(ierr)?; 485 Ok(Self { 486 qf_core: QFunctionCore { ceed, ptr }, 487 qf_ctx_ptr, 488 trampoline_data, 489 }) 490 } 491 492 /// Apply the action of a QFunction 493 /// 494 /// * `Q` - The number of quadrature points 495 /// * `input` - Array of input Vectors 496 /// * `output` - Array of output Vectors 497 /// 498 /// ``` 499 /// # use libceed::prelude::*; 500 /// # fn main() -> Result<(), libceed::CeedError> { 501 /// # let ceed = libceed::Ceed::default_init(); 502 /// let mut user_f = |[u, weights, ..]: QFunctionInputs, [v, ..]: QFunctionOutputs| { 503 /// // Iterate over quadrature points 504 /// v.iter_mut() 505 /// .zip(u.iter().zip(weights.iter())) 506 /// .for_each(|(v, (u, w))| *v = u * w); 507 /// 508 /// // Return clean error code 509 /// 0 510 /// }; 511 /// 512 /// let qf = ceed 513 /// .q_function_interior(1, Box::new(user_f))? 514 /// .input("u", 1, EvalMode::Interp)? 515 /// .input("weights", 1, EvalMode::Weight)? 516 /// .output("v", 1, EvalMode::Interp)?; 517 /// 518 /// const Q: usize = 8; 519 /// let mut w = [0.; Q]; 520 /// let mut u = [0.; Q]; 521 /// let mut v = [0.; Q]; 522 /// 523 /// for i in 0..Q { 524 /// let x = 2. * (i as Scalar) / ((Q as Scalar) - 1.) - 1.; 525 /// u[i] = 2. + 3. * x + 5. * x * x; 526 /// w[i] = 1. - x * x; 527 /// v[i] = u[i] * w[i]; 528 /// } 529 /// 530 /// let uu = ceed.vector_from_slice(&u)?; 531 /// let ww = ceed.vector_from_slice(&w)?; 532 /// let mut vv = ceed.vector(Q)?; 533 /// vv.set_value(0.0); 534 /// { 535 /// let input = vec![uu, ww]; 536 /// let mut output = vec![vv]; 537 /// qf.apply(Q, &input, &output)?; 538 /// vv = output.remove(0); 539 /// } 540 /// 541 /// vv.view() 542 /// .iter() 543 /// .zip(v.iter()) 544 /// .for_each(|(computed, actual)| { 545 /// assert_eq!( 546 /// *computed, *actual, 547 /// "Incorrect value in QFunction application" 548 /// ); 549 /// }); 550 /// # Ok(()) 551 /// # } 552 /// ``` 553 pub fn apply(&self, Q: usize, u: &[Vector], v: &[Vector]) -> crate::Result<i32> { 554 self.qf_core.apply(Q, u, v) 555 } 556 557 /// Add a QFunction input 558 /// 559 /// * `fieldname` - Name of QFunction field 560 /// * `size` - Size of QFunction field, `(ncomp * dim)` for `Grad` or 561 /// `(ncomp * 1)` for `None`, `Interp`, and `Weight` 562 /// * `emode` - `EvalMode::None` to use values directly, `EvalMode::Interp` 563 /// to use interpolated values, `EvalMode::Grad` to use 564 /// gradients, `EvalMode::Weight` to use quadrature weights 565 /// 566 /// ``` 567 /// # use libceed::prelude::*; 568 /// # fn main() -> Result<(), libceed::CeedError> { 569 /// # let ceed = libceed::Ceed::default_init(); 570 /// let mut user_f = |[u, weights, ..]: QFunctionInputs, [v, ..]: QFunctionOutputs| { 571 /// // Iterate over quadrature points 572 /// v.iter_mut() 573 /// .zip(u.iter().zip(weights.iter())) 574 /// .for_each(|(v, (u, w))| *v = u * w); 575 /// 576 /// // Return clean error code 577 /// 0 578 /// }; 579 /// 580 /// let mut qf = ceed.q_function_interior(1, Box::new(user_f))?; 581 /// 582 /// qf = qf.input("u", 1, EvalMode::Interp)?; 583 /// qf = qf.input("weights", 1, EvalMode::Weight)?; 584 /// # Ok(()) 585 /// # } 586 /// ``` 587 pub fn input( 588 mut self, 589 fieldname: &str, 590 size: usize, 591 emode: crate::EvalMode, 592 ) -> crate::Result<Self> { 593 let name_c = CString::new(fieldname).expect("CString::new failed"); 594 let idx = self.trampoline_data.number_inputs; 595 self.trampoline_data.input_sizes[idx] = size; 596 self.trampoline_data.number_inputs += 1; 597 let (size, emode) = ( 598 i32::try_from(size).unwrap(), 599 emode as bind_ceed::CeedEvalMode, 600 ); 601 let ierr = unsafe { 602 bind_ceed::CeedQFunctionAddInput(self.qf_core.ptr, name_c.as_ptr(), size, emode) 603 }; 604 self.qf_core.ceed.check_error(ierr)?; 605 Ok(self) 606 } 607 608 /// Add a QFunction output 609 /// 610 /// * `fieldname` - Name of QFunction field 611 /// * `size` - Size of QFunction field, `(ncomp * dim)` for `Grad` or 612 /// `(ncomp * 1)` for `None` and `Interp` 613 /// * `emode` - `EvalMode::None` to use values directly, `EvalMode::Interp` 614 /// to use interpolated values, `EvalMode::Grad` to use 615 /// gradients 616 /// 617 /// ``` 618 /// # use libceed::prelude::*; 619 /// # fn main() -> Result<(), libceed::CeedError> { 620 /// # let ceed = libceed::Ceed::default_init(); 621 /// let mut user_f = |[u, weights, ..]: QFunctionInputs, [v, ..]: QFunctionOutputs| { 622 /// // Iterate over quadrature points 623 /// v.iter_mut() 624 /// .zip(u.iter().zip(weights.iter())) 625 /// .for_each(|(v, (u, w))| *v = u * w); 626 /// 627 /// // Return clean error code 628 /// 0 629 /// }; 630 /// 631 /// let mut qf = ceed.q_function_interior(1, Box::new(user_f))?; 632 /// 633 /// qf.output("v", 1, EvalMode::Interp)?; 634 /// # Ok(()) 635 /// # } 636 /// ``` 637 pub fn output( 638 mut self, 639 fieldname: &str, 640 size: usize, 641 emode: crate::EvalMode, 642 ) -> crate::Result<Self> { 643 let name_c = CString::new(fieldname).expect("CString::new failed"); 644 let idx = self.trampoline_data.number_outputs; 645 self.trampoline_data.output_sizes[idx] = size; 646 self.trampoline_data.number_outputs += 1; 647 let (size, emode) = ( 648 i32::try_from(size).unwrap(), 649 emode as bind_ceed::CeedEvalMode, 650 ); 651 let ierr = unsafe { 652 bind_ceed::CeedQFunctionAddOutput(self.qf_core.ptr, name_c.as_ptr(), size, emode) 653 }; 654 self.qf_core.ceed.check_error(ierr)?; 655 Ok(self) 656 } 657 } 658 659 // ----------------------------------------------------------------------------- 660 // QFunction 661 // ----------------------------------------------------------------------------- 662 impl<'a> QFunctionByName<'a> { 663 // Constructor 664 pub fn create(ceed: &'a crate::Ceed, name: &str) -> crate::Result<Self> { 665 let name_c = CString::new(name).expect("CString::new failed"); 666 let mut ptr = std::ptr::null_mut(); 667 let ierr = unsafe { 668 bind_ceed::CeedQFunctionCreateInteriorByName(ceed.ptr, name_c.as_ptr(), &mut ptr) 669 }; 670 ceed.check_error(ierr)?; 671 Ok(Self { 672 qf_core: QFunctionCore { ceed, ptr }, 673 }) 674 } 675 676 /// Apply the action of a QFunction 677 /// 678 /// * `Q` - The number of quadrature points 679 /// * `input` - Array of input Vectors 680 /// * `output` - Array of output Vectors 681 /// 682 /// ``` 683 /// # use libceed::prelude::*; 684 /// # fn main() -> Result<(), libceed::CeedError> { 685 /// # let ceed = libceed::Ceed::default_init(); 686 /// const Q: usize = 8; 687 /// let qf_build = ceed.q_function_interior_by_name("Mass1DBuild")?; 688 /// let qf_mass = ceed.q_function_interior_by_name("MassApply")?; 689 /// 690 /// let mut j = [0.; Q]; 691 /// let mut w = [0.; Q]; 692 /// let mut u = [0.; Q]; 693 /// let mut v = [0.; Q]; 694 /// 695 /// for i in 0..Q { 696 /// let x = 2. * (i as Scalar) / ((Q as Scalar) - 1.) - 1.; 697 /// j[i] = 1.; 698 /// w[i] = 1. - x * x; 699 /// u[i] = 2. + 3. * x + 5. * x * x; 700 /// v[i] = w[i] * u[i]; 701 /// } 702 /// 703 /// let jj = ceed.vector_from_slice(&j)?; 704 /// let ww = ceed.vector_from_slice(&w)?; 705 /// let uu = ceed.vector_from_slice(&u)?; 706 /// let mut vv = ceed.vector(Q)?; 707 /// vv.set_value(0.0); 708 /// let mut qdata = ceed.vector(Q)?; 709 /// qdata.set_value(0.0); 710 /// 711 /// { 712 /// let mut input = vec![jj, ww]; 713 /// let mut output = vec![qdata]; 714 /// qf_build.apply(Q, &input, &output)?; 715 /// qdata = output.remove(0); 716 /// } 717 /// 718 /// { 719 /// let mut input = vec![qdata, uu]; 720 /// let mut output = vec![vv]; 721 /// qf_mass.apply(Q, &input, &output)?; 722 /// vv = output.remove(0); 723 /// } 724 /// 725 /// vv.view() 726 /// .iter() 727 /// .zip(v.iter()) 728 /// .for_each(|(computed, actual)| { 729 /// assert_eq!( 730 /// *computed, *actual, 731 /// "Incorrect value in QFunction application" 732 /// ); 733 /// }); 734 /// # Ok(()) 735 /// # } 736 /// ``` 737 pub fn apply(&self, Q: usize, u: &[Vector], v: &[Vector]) -> crate::Result<i32> { 738 self.qf_core.apply(Q, u, v) 739 } 740 } 741 742 // ----------------------------------------------------------------------------- 743