| 608 | } |
| 609 | |
| 610 | void geom_epsilon::eff_chi1inv_row(meep::component c, double chi1inv_row[3], const meep::volume &v, |
| 611 | double tol, int maxeval) { |
| 612 | const geometric_object *o; |
| 613 | material_type mat, mat_behind; |
| 614 | symmetric_matrix meps, meps_inv; |
| 615 | vector3 p, shiftby, normal; |
| 616 | bool destroy_material = false; |
| 617 | |
| 618 | if (maxeval == 0 || !get_front_object(v, geometry_tree, p, &o, shiftby, mat, mat_behind)) { |
| 619 | noavg: |
| 620 | destroy_material = get_material_pt(mat, v.center()); |
| 621 | trivial: |
| 622 | material_epsmu(meep::type(c), mat, &meps, &meps_inv); |
| 623 | switch (component_direction(c)) { |
| 624 | case meep::X: |
| 625 | case meep::R: |
| 626 | chi1inv_row[0] = meps_inv.m00; |
| 627 | chi1inv_row[1] = meps_inv.m01; |
| 628 | chi1inv_row[2] = meps_inv.m02; |
| 629 | break; |
| 630 | case meep::Y: |
| 631 | case meep::P: |
| 632 | chi1inv_row[0] = meps_inv.m01; |
| 633 | chi1inv_row[1] = meps_inv.m11; |
| 634 | chi1inv_row[2] = meps_inv.m12; |
| 635 | break; |
| 636 | case meep::Z: |
| 637 | chi1inv_row[0] = meps_inv.m02; |
| 638 | chi1inv_row[1] = meps_inv.m12; |
| 639 | chi1inv_row[2] = meps_inv.m22; |
| 640 | break; |
| 641 | case meep::NO_DIRECTION: chi1inv_row[0] = chi1inv_row[1] = chi1inv_row[2] = 0; |
| 642 | } |
| 643 | if (destroy_material) material_type_destroy(mat); |
| 644 | return; |
| 645 | } |
| 646 | |
| 647 | // FIXME: reimplement support for fallback integration, without |
| 648 | // messing up anisotropic support |
| 649 | // if (!get_front_object(v, geometry_tree, |
| 650 | // p, &o, shiftby, mat, mat_behind)) { |
| 651 | // fallback_chi1inv_row(c, chi1inv_row, v, tol, maxeval); |
| 652 | // return; |
| 653 | // } |
| 654 | |
| 655 | /* check for trivial case of only one object/material */ |
| 656 | if (material_type_equal(&mat, &mat_behind)) goto trivial; |
| 657 | |
| 658 | // it doesn't make sense to average metals (electric or magnetic) |
| 659 | if (is_metal(meep::type(c), &mat) || is_metal(meep::type(c), &mat_behind)) goto noavg; |
| 660 | |
| 661 | normal = unit_vector3(normal_to_fixed_object(vector3_minus(p, shiftby), *o)); |
| 662 | if (normal.x == 0 && normal.y == 0 && normal.z == 0) |
| 663 | goto noavg; // couldn't get normal vector for this point, punt |
| 664 | geom_box pixel = gv2box(v); |
| 665 | pixel.low = vector3_minus(pixel.low, shiftby); |
| 666 | pixel.high = vector3_minus(pixel.high, shiftby); |
| 667 |
nothing calls this directly
no test coverage detected