Add a test in PolynomialP

tests/test_PolynomialP.cpp

@@ -157,6 +157,7 @@ TEST_CASE("PolynomialP", "[analysis]")
     std::array<TinyVector<2>, 4> positions;
     std::array<TinyVector<2>, 3> positions2;
     std::array<TinyVector<2>, 2> positions4;
+    std::array<TinyVector<2>, 4> positions3;
     positions[0]  = TinyVector<2>{0, 0};
     positions[1]  = TinyVector<2>{0, 0.5};
     positions[2]  = TinyVector<2>{0.3, 0.7};
@@ -166,6 +167,10 @@ TEST_CASE("PolynomialP", "[analysis]")
     positions2[2] = TinyVector<2>{0.3, 0.7};
     positions4[0] = TinyVector<2>{0, 0.5};
     positions4[1] = TinyVector<2>{0.3, 0.7};
+    positions3[0] = TinyVector<2>{0, 0};
+    positions3[1] = TinyVector<2>{1, 0};
+    positions3[2] = TinyVector<2>{1, 1};
+    positions3[3] = TinyVector<2>{0, 1};
 
     PolynomialP<2, 2> P(coef);
     auto p1 = [](const TinyVector<2>& X) {
@@ -196,9 +201,15 @@ TEST_CASE("PolynomialP", "[analysis]")
     for (size_t i = 0; i < l1.numberOfPoints(); ++i) {
       value4 += l1.weight(i) * u.velocityNorm() * p1(u(l1.point(i)));
     }
+    SquareTransformation<2> s3{positions3[0], positions3[1], positions3[2], positions3[3]};
+    auto value3 = weight_list[0] * s3.jacobianDeterminant(point_list[0]) * p1(s3(point_list[0]));
+    for (size_t i = 1; i < weight_list.size(); ++i) {
+      value3 += weight_list[i] * s3.jacobianDeterminant(point_list[i]) * p1(s3(point_list[i]));
+    }
 
     REQUIRE(value == Catch::Approx(integrate(P, positions)));
     REQUIRE(value2 == Catch::Approx(integrate(P, positions2)));
+    REQUIRE(value3 == Catch::Approx(integrate(P, positions3)));
     REQUIRE(value4 == Catch::Approx(integrate(P, positions4)));
   }