New implementation of interface sources.

include/maxwell/maxwell_common.h

@@ -265,6 +265,47 @@ eval_interface_sources(const model& mod, const parameter_loader& mpl,
     }
 }
 
+/* Field must be either 'field' or 'pinned_field'. Too early to use a concept. */
+template<typename State, typename Field>
+void
+eval_interface_sources_new(const model& mod, const parameter_loader& mpl,
+    State& state, Field& srcfield)
+{
+    auto& is = mod.get_interfaces();
+    assert( srcfield.num_dofs == mod.num_fluxes() );
+
+    for (auto& [si, sd] : is)
+    {
+        const auto& e = mod.entity_at(sd.adjacent_volume_entity);
+        //auto etag = e.gmsh_tag();
+        auto mtag = sd.gmsh_tag();
+        switch ( mpl.interface_type( mtag ) )
+        {
+            case interface_condition::UNSPECIFIED:
+            case interface_condition::E_FIELD:
+            case interface_condition::H_FIELD:
+                break;
+
+            case interface_condition::SURFACE_CURRENT: {
+                for (auto& [local, global] : sd.face_numbers)
+                {
+                    vec3d n = state.eds[e.number()].normals.row(local);
+                    auto fJ = [&](const point_3d& pt) -> vec3d {
+                        vec3d J = mpl.eval_interface_source(mtag, pt, state.curr_time);
+                        return n.cross(J);
+                    };
+
+                    matxd prjJ = e.project_on_face(local, fJ);
+                    auto num_bf = prjJ.rows();
+                    srcfield.Hx.segment(global*num_bf, num_bf) = prjJ.col(0);
+                    srcfield.Hy.segment(global*num_bf, num_bf) = prjJ.col(1);
+                    srcfield.Hz.segment(global*num_bf, num_bf) = prjJ.col(2);
+                }
+            } break; /* case interface_condition::SURFACE_CURRENT */
+        }
+    }
+}
+
 #ifdef USE_MPI
 std::vector<double> compute_cell_ranks(const model&);
 #endif /* USE_MPI */

src/libgmshdg/gmsh_io.cpp

@@ -576,10 +576,10 @@ model::map_boundaries_new(void)
         sd.sort();
     }
 
-    //std::cout << "Interfaces:" << std::endl;
+    std::cout << "Interfaces:" << std::endl;
     for (auto& [si, sd] : interfaces)
     {
-        //std::cout << si << " " << sd.face_numbers.size() << std::endl;
+        std::cout << si << " " << sd.face_numbers.size() << std::endl;
         sd.sort();
     }
 }