Skip to content

Instantly share code, notes, and snippets.

@mdavezac

mdavezac/SoptShopTalk.ipynb

Last active October 26, 2016 10:11
Show Gist options
  • Save mdavezac/21b39888fe877aa39dd805984f961a89 to your computer and use it in GitHub Desktop.
Save mdavezac/21b39888fe877aa39dd805984f961a89 to your computer and use it in GitHub Desktop.
Download ZIP
Short talk describing the sopt (https://github.com/basp-group/sopt) compressed sensing library
Display the source blob
Display the rendered blob
Raw
SoptShopTalk.ipynb
{
"cells": [
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"# Sopt, a compressed sensing library\n",
"\n",
"- General algorithms for general convex minimization problems\n",
"- C++11\n",
"- Open source\n",
"- version controlled\n",
"- Unit-tested\n",
"- Continuously tested\n",
"- Now has a [shop-talk](https://gist.github.com/mdavezac/21b39888fe877aa39dd805984f961a89)! downloadable at http://bit.ly/2e9T8Md"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": false,
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"# The technical bits"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": false,
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"## Where to get it\n",
"\n",
"Available on github (https://github.com/basp-group/sopt)\n",
"under the GPL2:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": false,
"slideshow": {
"slide_type": "fragment"
}
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"Cloning into 'sopt'...\n"
]
}
],
"source": [
"%%bash\n",
"[ ! -e sopt ] && git clone https://github.com/basp-group/sopt"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": false,
"slideshow": {
"slide_type": "fragment"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"sopt:\r\n",
"CMakeLists.txt README.md \u001b[34mcpp\u001b[m\u001b[m/ \u001b[34mmatlab\u001b[m\u001b[m/\r\n",
"LICENSE.txt \u001b[34mcmake_files\u001b[m\u001b[m/ \u001b[34mimages\u001b[m\u001b[m/ \u001b[34mpython\u001b[m\u001b[m/\r\n",
"\r\n",
"sopt/cpp:\r\n",
"CMakeLists.txt \u001b[34mexamples\u001b[m\u001b[m/ \u001b[34msopt\u001b[m\u001b[m/ \u001b[34mtools_for_tests\u001b[m\u001b[m/\r\n",
"\u001b[34mbenchmarks\u001b[m\u001b[m/ \u001b[34mregressions\u001b[m\u001b[m/ \u001b[34mtests\u001b[m\u001b[m/\r\n"
]
}
],
"source": [
"ls sopt sopt/cpp"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"## installation and testing"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": false,
"slideshow": {
"slide_type": "subslide"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"-- The CXX compiler identification is AppleClang 8.0.0.8000038\n",
"-- Check for working CXX compiler: /Applications/Xcode.app/Contents/Developer/Toolchains/XcodeDefault.xctoolchain/usr/bin/clang++\n",
"-- Check for working CXX compiler: /Applications/Xcode.app/Contents/Developer/Toolchains/XcodeDefault.xctoolchain/usr/bin/clang++ -- works\n",
"-- Detecting CXX compiler ABI info\n",
"-- Detecting CXX compiler ABI info - done\n",
"-- Detecting CXX compile features\n",
"-- Detecting CXX compile features - done\n",
"-- [GreatCMakeCookOff] not found. Will attempt to clone it.\n",
"-- Found Git: /usr/local/bin/git (found version \"2.10.1\") \n",
"-- [GreatCMakeCookOff] downloaded to /Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/build/external/src/GreatCMakeCookOff\n",
"-- Performing Test has_std_cpp11\n",
"-- Performing Test has_std_cpp11 - Success\n",
"-- Performing Test has_std_cpp0x\n",
"-- Performing Test has_std_cpp0x - Success\n",
"-- Checking C++11 support for \"unique_ptr\"\n",
"-- Checking C++11 support for \"unique_ptr\" -- works\n",
"-- Checking C++11 support for \"nullptr\" (N2431)\n",
"-- Checking C++11 support for \"nullptr\" (N2431) -- works\n",
"-- Checking C++11 support for \"override\"\n",
"-- Checking C++11 support for \"override\" -- works\n",
"-- Checking C++11 support for \"constructor_delegate\"\n",
"-- Checking C++11 support for \"constructor_delegate\" -- works\n",
"-- Found Eigen3\n",
"-- Found spdlog\n",
"-- Found TIFF: /usr/local/lib/libtiff.dylib (found version \"4.0.6\") \n",
"-- Try OpenMP CXX flag = [-fopenmp=libomp]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Try OpenMP CXX flag = [ ]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Try OpenMP CXX flag = [-fopenmp]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Try OpenMP CXX flag = [/openmp]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Try OpenMP CXX flag = [-Qopenmp]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Try OpenMP CXX flag = [-openmp]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Try OpenMP CXX flag = [-xopenmp]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Try OpenMP CXX flag = [+Oopenmp]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Try OpenMP CXX flag = [-qsmp]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Try OpenMP CXX flag = [-mp]\n",
"-- Performing Test OpenMP_FLAG_DETECTED\n",
"-- Performing Test OpenMP_FLAG_DETECTED - Failed\n",
"-- Could NOT find OpenMP (missing: OpenMP_CXX_FLAGS) \n",
"-- Could not find OpenMP. Compiling without.\n",
"-- Found Catch: /Users/mdavezac/spack/opt/spack/darwin-sierra-x86_64/clang-8.0.0-apple/Catch-1.3.0-ktlhhdjh6f4nwlu5kopq2olx344r22ht/include \n",
"-- Configuring done\n",
"-- Generating done\n",
"-- Build files have been written to: /Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/build\n",
"Scanning dependencies of target lookup_dependencies\n",
"[ 0%] Built target lookup_dependencies\n",
"Scanning dependencies of target common_catch_main_object\n",
"[ 3%] Building CXX object cpp/tests/CMakeFiles/common_catch_main_object.dir/common_catch_main.cc.o\n",
"Scanning dependencies of target sopt\n",
"[ 12%] Building CXX object cpp/sopt/CMakeFiles/sopt.dir/wavelets/wavelets.cc.o\n",
"[ 12%] Building CXX object cpp/sopt/CMakeFiles/sopt.dir/utilities.cc.o\n",
"[ 12%] Building CXX object cpp/sopt/CMakeFiles/sopt.dir/wavelets/wavelet_data.cc.o\n",
"[ 16%] Linking CXX shared library libsopt.dylib\n",
"[ 16%] Built target sopt\n",
"[ 16%] Built target common_catch_main_object\n",
"Scanning dependencies of target test_conjugate_gradient\n",
"Scanning dependencies of target test_linear_transform\n",
"Scanning dependencies of target test_proximal\n",
"Scanning dependencies of target test_sara\n",
"[ 29%] Building CXX object cpp/tests/CMakeFiles/test_proximal.dir/proximal.cc.o\n",
"[ 29%] Building CXX object cpp/tests/CMakeFiles/test_conjugate_gradient.dir/conjugate_gradient.cc.o\n",
"[ 29%] Building CXX object cpp/tests/CMakeFiles/test_sara.dir/sara.cc.o\n",
"[ 29%] Building CXX object cpp/tests/CMakeFiles/test_linear_transform.dir/linear_transform.cc.o\n",
"Scanning dependencies of target test_sdmm_warm_start\n",
"Scanning dependencies of target test_padmm_warm_start\n",
"Scanning dependencies of target test_padmm\n",
"Scanning dependencies of target test_reweighted\n",
"[ 35%] Building CXX object cpp/tests/CMakeFiles/test_padmm.dir/padmm.cc.o\n",
"[ 35%] Building CXX object cpp/tests/CMakeFiles/test_sdmm_warm_start.dir/sdmm_warm_start.cc.o\n",
"[ 38%] Building CXX object cpp/tests/CMakeFiles/test_padmm_warm_start.dir/padmm_warm_start.cc.o\n",
"[ 41%] Building CXX object cpp/tests/CMakeFiles/test_reweighted.dir/reweighted.cc.o\n",
"[ 45%] Linking CXX executable test_reweighted\n",
"[ 45%] Built target test_reweighted\n",
"Scanning dependencies of target test_maths\n",
"[ 48%] Building CXX object cpp/tests/CMakeFiles/test_maths.dir/maths.cc.o\n",
"[ 51%] Linking CXX executable test_linear_transform\n",
"[ 51%] Built target test_linear_transform\n",
"Scanning dependencies of target test_wavelets\n",
"[ 54%] Building CXX object cpp/tests/CMakeFiles/test_wavelets.dir/wavelets.cc.o\n",
"[ 58%] Linking CXX executable test_padmm_warm_start\n",
"[ 58%] Built target test_padmm_warm_start\n",
"Scanning dependencies of target test_power_method\n",
"[ 61%] Building CXX object cpp/tests/CMakeFiles/test_power_method.dir/power_method.cc.o\n",
"[ 64%] Linking CXX executable test_conjugate_gradient\n",
"[ 64%] Built target test_conjugate_gradient\n",
"[ 67%] Linking CXX executable test_sdmm_warm_start\n",
"Scanning dependencies of target test_sdmm\n",
"[ 70%] Building CXX object cpp/tests/CMakeFiles/test_sdmm.dir/sdmm.cc.o\n",
"[ 70%] Built target test_sdmm_warm_start\n",
"Scanning dependencies of target test_wrapper\n",
"[ 74%] Building CXX object cpp/tests/CMakeFiles/test_wrapper.dir/wrapper.cc.o\n",
"[ 77%] Linking CXX executable test_padmm\n",
"[ 77%] Built target test_padmm\n",
"[ 80%] Linking CXX executable test_wrapper\n",
"[ 83%] Linking CXX executable test_proximal\n",
"[ 83%] Built target test_wrapper\n",
"[ 83%] Built target test_proximal\n",
"[ 87%] Linking CXX executable test_maths\n",
"[ 87%] Built target test_maths\n",
"[ 90%] Linking CXX executable test_sara\n",
"[ 90%] Built target test_sara\n",
"[ 93%] Linking CXX executable test_wavelets\n",
"[ 96%] Linking CXX executable test_power_method\n",
"[ 96%] Built target test_wavelets\n",
"[ 96%] Built target test_power_method\n",
"[100%] Linking CXX executable test_sdmm\n",
"[100%] Built target test_sdmm\n"
]
}
],
"source": [
"%%bash\n",
"mkdir -p sopt/build && cd sopt/build\n",
"cmake -DCMAKE_BUILD_TYPE=Release -Dexamples=OFF ..\n",
"make -j 8"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": false,
"slideshow": {
"slide_type": "subslide"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Running tests...\n",
"Test project /Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/build\n",
" Start 1: test_wavelets\n",
" 1/14 Test #1: test_wavelets .................... Passed 0.28 sec\n",
" Start 2: test_sara\n",
" 2/14 Test #2: test_sara ........................ Passed 0.05 sec\n",
" Start 3: test_maths\n",
" 3/14 Test #3: test_maths ....................... Passed 0.01 sec\n",
" Start 4: test_wrapper\n",
" 4/14 Test #4: test_wrapper ..................... Passed 0.00 sec\n",
" Start 5: test_conjugate_gradient\n",
" 5/14 Test #5: test_conjugate_gradient .......... Passed 0.00 sec\n",
" Start 6: test_linear_transform\n",
" 6/14 Test #6: test_linear_transform ............ Passed 0.00 sec\n",
" Start 7: test_sdmm\n",
" 7/14 Test #7: test_sdmm ........................ Passed 0.08 sec\n",
" Start 8: test_sdmm_warm_start\n",
" 8/14 Test #8: test_sdmm_warm_start ............. Passed 0.01 sec\n",
" Start 9: test_proximal\n",
" 9/14 Test #9: test_proximal .................... Passed 0.01 sec\n",
" Start 10: seeded_proximal\n",
"10/14 Test #10: seeded_proximal .................. Passed 0.01 sec\n",
" Start 11: test_padmm\n",
"11/14 Test #11: test_padmm ....................... Passed 0.01 sec\n",
" Start 12: test_padmm_warm_start\n",
"12/14 Test #12: test_padmm_warm_start ............ Passed 0.01 sec\n",
" Start 13: test_reweighted\n",
"13/14 Test #13: test_reweighted .................. Passed 0.01 sec\n",
" Start 14: test_power_method\n",
"14/14 Test #14: test_power_method ................ Passed 0.00 sec\n",
"\n",
"100% tests passed, 0 tests failed out of 14\n",
"\n",
"Label Time Summary:\n",
"catch = 0.47 sec (14 tests)\n",
"\n",
"Total Test time (real) = 0.48 sec\n"
]
}
],
"source": [
"%%bash\n",
"cd sopt/build && make test"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"## How we know it works\n",
"\n",
"- Each functional bit of code is unit-tested"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false,
"slideshow": {
"slide_type": "fragment"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"CMakeLists.txt padmm.cc sara.cc\r\n",
"common_catch_main.cc padmm_warm_start.cc sdmm.cc\r\n",
"conjugate_gradient.cc power_method.cc sdmm_warm_start.cc\r\n",
"linear_transform.cc proximal.cc wavelets.cc\r\n",
"maths.cc reweighted.cc wrapper.cc\r\n"
]
}
],
"source": [
"ls sopt/cpp/tests"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"```cpp\n",
"TEST_CASE(\"Projector on positive quadrant\", \"[utility][project]\") {\n",
" using namespace sopt;\n",
"\n",
" SECTION(\"Real matrix\") {\n",
" Image<> input = Image<>::Ones(5, 5) + Image<>::Random(5, 5) * 0.55;\n",
" input(1, 1) *= -1;\n",
" input(3, 2) *= -1;\n",
"\n",
" auto const expr = positive_quadrant(input);\n",
" CHECK(expr(1, 1) == Approx(0));\n",
" CHECK(expr(3, 2) == Approx(0));\n",
" }\n",
"}\n",
"```"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"collapsed": false,
"slideshow": {
"slide_type": "fragment"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n",
"test_maths is a Catch v1.3.0 host application.\n",
"Run with -? for options\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Projector on positive quadrant\n",
" Real matrix\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:11\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:22: \n",
"PASSED:\n",
" CHECK( expr(1, 1) == Approx(0) )\n",
"with expansion:\n",
" 0.0 == Approx( 0.0 )\n",
"with messages:\n",
" input := 0.450009 0.690855 0.871852 1.03267 0.909235\n",
" 0.594692 -0.501749 1.02136 1.18826 1.20545\n",
" 1.28117 1.19675 1.36406 0.458468 1.09787\n",
" 0.954515 1.19723 -0.488029 0.871757 1.47348\n",
" 1.03604 1.47816 0.508808 0.523526 1.38078\n",
" expr := 0.450009 0.690855 0.871852 1.03267 0.909235\n",
" 0.594692 0 1.02136 1.18826 1.20545\n",
" 1.28117 1.19675 1.36406 0.458468 1.09787\n",
" 0.954515 1.19723 0 0.871757 1.47348\n",
" 1.03604 1.47816 0.508808 0.523526 1.38078\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:23: \n",
"PASSED:\n",
" CHECK( expr(3, 2) == Approx(0) )\n",
"with expansion:\n",
" 0.0 == Approx( 0.0 )\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:26: \n",
"PASSED:\n",
" CHECK( value(1, 1) == Approx(0) )\n",
"with expansion:\n",
" 0.0 == Approx( 0.0 )\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:27: \n",
"PASSED:\n",
" CHECK( value(3, 2) == Approx(0) )\n",
"with expansion:\n",
" 0.0 == Approx( 0.0 )\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:30: \n",
"PASSED:\n",
" CHECK( value.isApprox(input) )\n",
"with expansion:\n",
" true\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Projector on positive quadrant\n",
" Complex matrix\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:11\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:41: \n",
"PASSED:\n",
" CHECK( expr.imag().isApprox(Image<>::Zero(5, 5)) )\n",
"with expansion:\n",
" true\n",
"with messages:\n",
" input := (1.02962,-0.448839) (0.811058,0.145902) (0.529954,0.144798)\n",
" (0.633158,-0.0148309) (0.993374,-0.257241)\n",
" (1.16931,-0.0924007) (-1.28205,0.540141) (1.42318,-0.250019) (1\n",
" .43742,0.450129) (0.549806,0.492541)\n",
" (1.22131,0.451353) (0.851873,-0.278257) (0.930053,0.293144) (0\n",
" .516621,0.445118) (0.531124,0.000777804)\n",
" (1.28842,-0.261302) (1.53081,0.244926) (-0.975505,-0.288448) (1\n",
" .00498,0.0179212) (0.872556,-0.24521)\n",
" (0.502211,0.25969) (1.27869,0.16667) (0.752398,-0.154809) (0\n",
" .800936,0.535306) (1.4552,0.0327221)\n",
" expr := (1.02962,0) (0.811058,0) (0.529954,0) (0.633158,0) (0.993374,0)\n",
" (1.16931,0) (0,0) (1.42318,0) (1.43742,0) (0.549806,0)\n",
" (1.22131,0) (0.851873,0) (0.930053,0) (0.516621,0) (0.531124,0)\n",
" (1.28842,0) (1.53081,0) (0,0) (1.00498,0) (0.872556,0)\n",
" (0.502211,0) (1.27869,0) (0.752398,0) (0.800936,0) (1.4552,0)\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:44: \n",
"PASSED:\n",
" CHECK( value.real()(1, 1) == Approx(0) )\n",
"with expansion:\n",
" 0.0 == Approx( 0.0 )\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:45: \n",
"PASSED:\n",
" CHECK( value.real()(3, 2) == Approx(0) )\n",
"with expansion:\n",
" 0.0 == Approx( 0.0 )\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:48: \n",
"PASSED:\n",
" CHECK( value.real().isApprox(input.real()) )\n",
"with expansion:\n",
" true\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:49: \n",
"PASSED:\n",
" CHECK( value.imag().isApprox(0e0 * input.real()) )\n",
"with expansion:\n",
" true\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Weighted l1 norm\n",
" Real valued\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:53\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:60: \n",
"PASSED:\n",
" CHECK( sopt::l1_norm(input, weight) == Approx(5 + 12 + 21 + 32) )\n",
"with expansion:\n",
" 70.0 == Approx( 70.0 )\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Weighted l1 norm\n",
" Complex valued\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:53\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:66: \n",
"PASSED:\n",
" CHECK( sopt::l1_norm(input, weight) == Approx(std::sqrt(2) * (5 + 12 + 21 + 32)) )\n",
"with expansion:\n",
" 98.9949493661 == Approx( 98.9949493661 )\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Soft threshhold\n",
" Single-valued threshhold\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:70\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:76: \n",
"PASSED:\n",
" CHECK( sopt::soft_threshhold(input, 1.1e1)(0) == Approx(0) )\n",
"with expansion:\n",
" 0.0 == Approx( 0.0 )\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:77: \n",
"PASSED:\n",
" CHECK( not(sopt::soft_threshhold(input, 1.1e1)(1) == Approx(0)) )\n",
"with expansion:\n",
" true\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:85: \n",
"PASSED:\n",
" CHECK( (b - a) == Approx(2 * (c - b)) )\n",
"with expansion:\n",
" 4.5 == Approx( 4.5 )\n",
"with messages:\n",
" b - a := 4.5\n",
" c - b := 2.25\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Soft threshhold\n",
" Multi-values threshhold\n",
" Real input\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:70\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:96: \n",
"PASSED:\n",
" CHECK( actual(0) == 0e0 )\n",
"with expansion:\n",
" 0.0 == 0.0\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:97: \n",
"PASSED:\n",
" CHECK( actual(1) == input(1) - threshhold(1) )\n",
"with expansion:\n",
" 9.0 == 9.0\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:98: \n",
"PASSED:\n",
" CHECK( actual(2) == input(2) + threshhold(2) )\n",
"with expansion:\n",
" -29.0 == -29.0\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:99: \n",
"PASSED:\n",
" CHECK( actual(3) == input(3) - threshhold(3) )\n",
"with expansion:\n",
" 35.5 == 35.5\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:101: \n",
"PASSED:\n",
" CHECK_THROWS_AS( soft_threshhold(input, threshhold.head(2)) )\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Soft threshhold\n",
" Multi-values threshhold\n",
" Complex input\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:70\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:105: \n",
"PASSED:\n",
" CHECK( actual(0) == 0e0 )\n",
"with expansion:\n",
" (0,0) == 0.0\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:106: \n",
"PASSED:\n",
" CHECK( actual(1) == input(1) - threshhold(1) )\n",
"with expansion:\n",
" (9,0) == 9.0\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:107: \n",
"PASSED:\n",
" CHECK( actual(2) == input(2) + threshhold(2) )\n",
"with expansion:\n",
" (-29,0) == -29.0\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:108: \n",
"PASSED:\n",
" CHECK( actual(3) == input(3) - threshhold(3) )\n",
"with expansion:\n",
" (35.5,0) == 35.5\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:110: \n",
"PASSED:\n",
" CHECK_THROWS_AS( soft_threshhold(input, threshhold.head(2)) )\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Sampling\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:115\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:123: \n",
"PASSED:\n",
" CHECK( down.size() == 4 )\n",
"with expansion:\n",
" 4 == 4\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:124: \n",
"PASSED:\n",
" CHECK( down(0) == input(1) )\n",
"with expansion:\n",
" 946997239 (0x387207f7)\n",
" ==\n",
" 946997239 (0x387207f7)\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:125: \n",
"PASSED:\n",
" CHECK( down(1) == input(3) )\n",
"with expansion:\n",
" 561597601 (0x21794ca1)\n",
" ==\n",
" 561597601 (0x21794ca1)\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:126: \n",
"PASSED:\n",
" CHECK( down(2) == input(6) )\n",
"with expansion:\n",
" -804521730 == -804521730\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:127: \n",
"PASSED:\n",
" CHECK( down(3) == input(5) )\n",
"with expansion:\n",
" 703982291 (0x29f5ead3)\n",
" ==\n",
" 703982291 (0x29f5ead3)\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:131: \n",
"PASSED:\n",
" CHECK( up(1) == input(1) )\n",
"with expansion:\n",
" 946997239 (0x387207f7)\n",
" ==\n",
" 946997239 (0x387207f7)\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:132: \n",
"PASSED:\n",
" CHECK( up(3) == input(3) )\n",
"with expansion:\n",
" 561597601 (0x21794ca1)\n",
" ==\n",
" 561597601 (0x21794ca1)\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:133: \n",
"PASSED:\n",
" CHECK( up(6) == input(6) )\n",
"with expansion:\n",
" -804521730 == -804521730\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:134: \n",
"PASSED:\n",
" CHECK( up(5) == input(5) )\n",
"with expansion:\n",
" 703982291 (0x29f5ead3)\n",
" ==\n",
" 703982291 (0x29f5ead3)\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:139: \n",
"PASSED:\n",
" CHECK( up == t_Vector::Zero(up.size()) )\n",
"with expansion:\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" ==\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
" 0\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Relative variation\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:142\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:146: \n",
"PASSED:\n",
" CHECK( not relvar(input) )\n",
"with expansion:\n",
" true\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:147: \n",
"PASSED:\n",
" CHECK( relvar(input) )\n",
"with expansion:\n",
" true\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:148: \n",
"PASSED:\n",
" CHECK( relvar(input + relvar.epsilon() * 0.5 / 6. * sopt::Array<>::Random(6)) )\n",
"with expansion:\n",
" true\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:149: \n",
"PASSED:\n",
" CHECK( not relvar(input + relvar.epsilon() * 1.1 * sopt::Array<>::Ones(6)) )\n",
"with expansion:\n",
" true\n",
"\n",
"-------------------------------------------------------------------------------\n",
"Standard deviation\n",
"-------------------------------------------------------------------------------\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:152\n",
"...............................................................................\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:160: \n",
"PASSED:\n",
" CHECK( std::abs(sopt::standard_deviation(input) - stddev) < 1e-8 )\n",
"with expansion:\n",
" 0.0 < 0.00000001\n",
"\n",
"/Users/mdavezac/workspaces/bico/src/soptshoptalk/sopt/cpp/tests/maths.cc:161: \n",
"PASSED:\n",
" CHECK( std::abs(sopt::standard_deviation(input.matrix()) - stddev) < 1e-8 )\n",
"with expansion:\n",
" 0.0 < 0.00000001\n",
"\n",
"===============================================================================\n",
"All tests passed (41 assertions in 6 test cases)\n",
"\n"
]
}
],
"source": [
"%%bash\n",
"sopt/build/cpp/tests/test_maths -s"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"## How we know it works\n",
"\n",
"- Each functional bit of code is unit-tested\n",
"- Continuously tested:\n",
"\n",
"![](jenkins.png)"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## Learning how to use it\n",
"\n",
"- take a look at the tests!\n",
"- take a look at this [talk](https://gist.github.com/mdavezac/21b39888fe877aa39dd805984f961a89)!\n",
"- examples! wow!"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"collapsed": false,
"slideshow": {
"slide_type": "fragment"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"CMakeLists.txt \u001b[34mproximal_admm\u001b[m\u001b[m/\r\n",
"conjugate_gradient.cc sara.cc\r\n",
"l1_norm.cc \u001b[34msdmm\u001b[m\u001b[m/\r\n",
"l1_proximal.cc soft_threshhold.cc\r\n",
"positive_quadrant_projection.cc wavelets.cc\r\n"
]
}
],
"source": [
"ls sopt/cpp/examples"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"```cpp\n",
"#include <sopt/maths.h>\n",
"#include <sopt/types.h>\n",
"\n",
"int main(int, char const **) {\n",
" sopt::Matrix<int> input(2, 2), weights(2, 2);\n",
" input << 1, -2, 3, -4;\n",
" weights << 5, 6, 7, 8;\n",
"\n",
" if(sopt::l1_norm(input, weights) != 1 * 5 + 2 * 6 + 3 * 7 + 4 * 8)\n",
" throw std::exception();\n",
"\n",
" return 0;\n",
"}\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"# Proximal Alternate Direction Method of Multipliers\n",
"\n",
"Solve $$\\min_{x\\in\\mathbb{R}^N, y\\in\\mathbb{R}^M} f(x) + g(z)\\quad\\text{s.t.}\\quad \\Phi x + z = y$$\n",
"\n",
"- $\\mathbf{L}$ is a linear transform/Measurement operator\n",
"- $y\\in\\mathbb{R}^M$ is a vector/target measurements\n",
"\n",
"The algorithm defines the *augmented Lagrangian* $\\mathcal{L}_\\gamma$\n",
"\n",
"$$\\mathcal{L}_\\gamma(x, y, \\lambda) = f(x) + g(z) + \\frac{\\lambda^H}{\\gamma}(\\Phi x + z - y) + \\frac{1}{2\\gamma}||\\Phi x + z - y||_2^2$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"### Algorithm\n",
"Loop until converged:\n",
"$$z^{(n+1)} = \\mathrm{prox}_{\\gamma g}(y - \\Phi x^{(n)} - \\lambda^{(n)})$$\n",
"$$x^{(n+1)} = \\mathrm{prox}_{\\frac{\\gamma}{\\nu} f}(x^{(n)}-\\mu\\Phi^Hx^{(n)} - y - z^{(n+1)})$$\n",
"$$\\lambda^{(n+1)} = \\lambda^{(n)} + \\beta(\\Phi x^{(n+1)} - y - z^{(n + 1)})$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"### Algorithm\n",
"Loop until converged:\n",
"$$z^{(n+1)} = \\mathrm{prox}_{\\gamma g}(y - \\Phi x^{(n)} - \\lambda^{(n)})$$\n",
"$$x^{(n+1)} = \\mathrm{prox}_{\\frac{\\gamma}{\\nu} f}(x^{(n)}-\\mu\\Phi^Hx^{(n)} - y - z^{(n+1)})$$\n",
"$$\\lambda^{(n+1)} = \\lambda^{(n)} + \\beta(\\Phi x^{(n+1)} - y - z^{(n + 1)})$$\n",
"\n",
"### Implementation\n",
"\n",
"```cpp\n",
"g_proximal(z, gamma(), -lambda - residual);\n",
"f_proximal(out, gamma() / nu(), out - Phi().adjoint() * (residual + lambda + z) / nu());\n",
"residual = Phi() * out - target();\n",
"lambda += lagrange_update_scale() * (residual + z);\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": true,
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"## Usage\n",
"\n",
"Lets solve $\\min_{x} ||x||_2 + ||x - x_1||_2$.\n",
"\n",
"The solution is any point on the segment $[0, x_1]$.\n"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"We will need a few headers:\n",
" \n",
"```cpp\n",
"#include <sopt/paddm.h>\n",
"#include <sopt/proximal.h>\n",
"#include <sopt/relative_variation.h>\n",
"#include <iostream>\n",
"```\n"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"Then we can define the proximals of `f` and `g`.\n",
"\n",
"Any function or lambda or functor with the right signature will work:\n",
"\n",
"```cpp\n",
"auto const prox_f = [](Vector<> &out, double gamma, Vector<> const &input) -> void {\n",
" auto const norm = x.stableNorm();\n",
" if(norm > t)\n",
" out = (1 - gamma / norm) * x;\n",
" else\n",
" out.fill(0); \n",
"};\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"Some helper functions are already defined in `sopt/proximal.h`:\n",
"\n",
"\n",
"```cpp\n",
"auto const N = 3;\n",
"auto const x0 = sopt::Vector<>::Random(N);\n",
"auto const prox_g = sopt::proximal::translate(sopt::proximal::EuclidianNorm(), -x0);\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"Then we defined the sensing operator $\\Phi$ and the target measurements $y$:\n",
" \n",
"```cpp\n",
"Vector<> const y = Vector<>::Zero(N);\n",
"Matrix<> const Phi = -Matrix<>::Identity(N, N);\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"We are now in a position to define the solver:\n",
"\n",
"```cpp\n",
"auto padmm = sopt::algorithm::ProximalADMM<double>(prox_f, prox_g, y)\n",
" .itermax(5000)\n",
" .is_converged(sopt::RelativeVariation<double>(1e-12))\n",
" .gamma(0.01)\n",
" .Phi(Phi);\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"And to use it:\n",
" \n",
"```cpp\n",
"auto const result = padmm();\n",
"std::cout << result.x.transpose() << std::endl;\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"`result` is a structure with the solution and information about the run:\n",
" \n",
"```cpp\n",
"struct {\n",
" //! Number of iterations\n",
" sopt::t_uint niters;\n",
" //! Wether convergence was achieved\n",
" bool good;\n",
" //! the residual from the last iteration\n",
" sopt::Vector<> residual;\n",
" //! The solution!\n",
" sopt::Vector<> x;\n",
"}\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"# Hands-on"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"## Putting the usage example in project of its own\n",
"\n",
"1. Create a CMake project that can find sopt and link to it:\n",
" In a new directory, create the file `CMakeLists.txt`\n",
" \n",
" ```CMake\n",
" cmake_minimum_required(VERSION 2.8)\n",
" project(Handson CXX)\n",
" \n",
" find_package(Sopt REQUIRED)\n",
" add_executable(main main.cc)\n",
" target_link_library(main sopt)\n",
" ```\n",
" \n",
"1. Put the example into `main.cc`:\n",
"\n",
" ```cpp\n",
" // headers\n",
" \n",
" int main(int, char**) {\n",
" // example\n",
" \n",
" return 0;\n",
" }\n",
" ```\n",
" \n",
"1. Build, compile and run:\n",
"\n",
" ```Bash\n",
" mkdir build\n",
" cd build\n",
" cmake ..\n",
" make\n",
" ./main\n",
" ```"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"## Hands-on 2\n",
"\n",
"Craft yourself something new to minimize: $\\min_x f(x) + g(x)$, $f(x) = ?$, $g(x) = ?$, including this proximal\n",
"\n",
"$$h(x) = x^T\\mathcal{A}x + b^Tx + c,\\quad A\\geq 0\\quad\\mapsto\\quad\n",
"\\mathrm{prox}_{\\gamma h}(x) = (I + tA)^{-1}(x - tb)$$\n",
"\n",
"And solve it using Sopt's PADMM method.\n",
"\n",
"----\n",
"Hint: Define new functions `prox_f` and `prox_g`, use them as input to `paddmm`. "
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"## Useful and less useful extras\n",
"\n",
"1. The verbosity of the logging/output is controlled by:\n",
" \n",
" ```cpp\n",
" // initialize logging in main\n",
" sopt::logging::initialize();\n",
" // control verbosity: trace, debug, info, warn, critical, off\n",
" sopt::logging::set_level(\"warn\");\n",
" ```"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"2. Linear transforms and their adjoint can be specified with a matrix or two functions (mapping a vector to another)\n",
"\n",
" ```cpp\n",
" auto const direct = [](Vector<> &out, Vector<> const &input) {\n",
" out = 2 * input;\n",
" };\n",
" auto const adjoint = [](Vector<> &out, Vector<> const &input) {\n",
" out = 1./2. * input;\n",
" };\n",
" auto const Phi = linear_transform(direct, adjoint, {1, 1, 0});\n",
" ```\n",
" \n",
" The last three numbers `(a, b, c)` help determine the size of the output vector from the size of the input vector: `size(output) = size(input) * a / b + c`.\n",
" If the size of the output `Nout` is fixed, then `{0, 0, Nout}` will do.\n",
" \n",
" The transform can be called on a `Vector<> x` as `Phi * x` and the adjoint as `Phi.adjoint() * x`.\n",
" "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"1. Some known proximal and known proximal transforms (translation, etc) are available in `cpp/sopt/proximal.h`\n",
"1. The proximal of the $x\\mapsto ||\\Phi x||_1$ is available in `sopt/cpp/l1_norm.h` "
]
}
],
"metadata": {
"celltoolbar": "Slideshow",
"kernelspec": {
"display_name": "Python 2",
"language": "python",
"name": "python2"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.12"
}
},
"nbformat": 4,
"nbformat_minor": 1
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment