diff --git a/unittest/force-styles/pair_style.cpp b/unittest/force-styles/pair_style.cpp
index e8da21e9fb..90bd646fc7 100644
--- a/unittest/force-styles/pair_style.cpp
+++ b/unittest/force-styles/pair_style.cpp
@@ -581,6 +581,51 @@ TEST(PairStyle, plain) {
     if (print_stats)
         std::cerr << "data_energy stats:" << stats << std::endl;
 
+    if (pair->respa_enable) {
+        if (!verbose) ::testing::internal::CaptureStdout();
+        cleanup_lammps(lmp,test_config);
+        lmp = init_lammps(argc,argv,test_config,false);
+        lmp->input->one("run_style respa 2 1 inner 1 4.0 5.0 outer 2");
+        run_lammps(lmp);
+        if (!verbose) ::testing::internal::GetCapturedStdout();
+
+        // need to relax error by a large amount with tabulation, since
+        // coul/long styles do not use tabulation in compute_inner()
+        // and compute_middle() so we get a significant deviation.
+        if (pair->ncoultablebits) epsilon *= 1.0e6;
+
+        f = lmp->atom->f;
+        tag=lmp->atom->tag;
+        pair = lmp->force->pair;
+        stats.reset();
+        for (int i=0; i < nlocal; ++i) {
+            EXPECT_FP_LE_WITH_EPS(f[i][0], f_run[tag[i]].x, 5*epsilon);
+            EXPECT_FP_LE_WITH_EPS(f[i][1], f_run[tag[i]].y, 5*epsilon);
+            EXPECT_FP_LE_WITH_EPS(f[i][2], f_run[tag[i]].z, 5*epsilon);
+        }
+        if (print_stats)
+            std::cerr << "run_forces  stats, r-RESPA:" << stats << std::endl;
+
+        stress = pair->virial;
+        stats.reset();
+        EXPECT_FP_LE_WITH_EPS(stress[0], test_config.run_stress.xx, epsilon);
+        EXPECT_FP_LE_WITH_EPS(stress[1], test_config.run_stress.yy, epsilon);
+        EXPECT_FP_LE_WITH_EPS(stress[2], test_config.run_stress.zz, epsilon);
+        EXPECT_FP_LE_WITH_EPS(stress[3], test_config.run_stress.xy, epsilon);
+        EXPECT_FP_LE_WITH_EPS(stress[4], test_config.run_stress.xz, epsilon);
+        EXPECT_FP_LE_WITH_EPS(stress[5], test_config.run_stress.yz, epsilon);
+        if (print_stats)
+            std::cerr << "run_stress  stats, r-RESPA:" << stats << std::endl;
+
+        stats.reset();
+        id = lmp->modify->find_compute("sum");
+        energy = lmp->modify->compute[id]->compute_scalar();
+        EXPECT_FP_LE_WITH_EPS(pair->eng_vdwl, test_config.run_vdwl, epsilon);
+        EXPECT_FP_LE_WITH_EPS(pair->eng_coul, test_config.run_coul, epsilon);
+        EXPECT_FP_LE_WITH_EPS((pair->eng_vdwl+pair->eng_coul),energy, epsilon);
+        if (print_stats)
+            std::cerr << "run_energy  stats, r-RESPA:" << stats << std::endl;
+    }
     if (!verbose) ::testing::internal::CaptureStdout();
     cleanup_lammps(lmp,test_config);
     if (!verbose) ::testing::internal::GetCapturedStdout();