diff --git a/fastlib2/contrib/nvasil/convex_nmf/build.py b/fastlib2/contrib/nvasil/convex_nmf/build.py
new file mode 100644
index 0000000000..60c3fbba4f
--- /dev/null
+++ b/fastlib2/contrib/nvasil/convex_nmf/build.py
@@ -0,0 +1,6 @@
+librule(name="nmflib",
+		headers=["nmf_objectives.h", "nmf_objectives_impl.h"],
+		deplibs=["fastlib:fastlib"] )
+binrule(name="nmf",
+		sources=["main.cc"],
+		deplibs=["fastlib:fastlib", ":nmflib", "contrib/nvasil/l_bfgs:l_bfgs"])
diff --git a/fastlib2/contrib/nvasil/convex_nmf/main.cc b/fastlib2/contrib/nvasil/convex_nmf/main.cc
new file mode 100644
index 0000000000..408331c785
--- /dev/null
+++ b/fastlib2/contrib/nvasil/convex_nmf/main.cc
@@ -0,0 +1,17 @@
+#include "fastlib/fastlib.h"
+#include "nmf_engine.h"
+
+int main(int argc, char *argv[]) {
+  fx_module *fx_root;
+	fx_root=fx_init(argc, argv, NULL);
+	fx_module nmf_module=fx_submodule(fx_root, "/engine");
+	NmfEngine engine;
+	engine.Init(nmf_mofule);
+  engine.ComputeNmf();
+	Matrix w_mat;
+	Matrix h_mat;
+  engine.GetW(&w_mat);
+	engine.GetH(&h_mat);
+  data::Save("W.csv", w_mat);
+	data::Save("H.csv", h_mat);
+}
diff --git a/fastlib2/contrib/nvasil/convex_nmf/nmf_engine.h b/fastlib2/contrib/nvasil/convex_nmf/nmf_engine.h
new file mode 100644
index 0000000000..31224a6d58
--- /dev/null
+++ b/fastlib2/contrib/nvasil/convex_nmf/nmf_engine.h
@@ -0,0 +1,80 @@
+#ifndef NMF_ENGINE_H_
+#define NMF_ENGINE_H_
+#include "fastlib/fastlib.h"
+#include "../l_bfgs/l_bfgs.h"
+#include "nmf_objectives.h"
+
+class NmfEngine {
+ public:
+	 void  Init(fx_module *module) {
+ 		 module_=module;
+	   std::string data_file=fx_param_str_req(module_, "data_file");
+     sdp_rank_=fx_param_int(module_,  "sdp_rank", 5);
+		 new_dim_=fx_param_int(module_, "new_dim", 3);
+		 Matrix data_mat;
+		 data::Load(data_file.c_str(), &data_mat);
+		 PreprcessData(data_mat);
+     opt_function_.Init(fx_submodule(module,"optfun"), 
+				 rows_, columns_, values_);
+     engine.Init(&opt_function, l_bfgs_node);
+	 }
+	 void Destruct() {
+	   
+	 };
+	 void ComputeNmf() {
+		  Matrix init_data;
+			opt_function_.GiveInitMatrix(&init_data);
+		  engine.set_coordinates(init_data);
+	    engine_.ComputeLocalOptimumBFGS();
+      Matrix result;
+      engine_.GetResults(&result);
+			w_mat_.Init(num_of_rows_, new_dim_);
+			h_mat_.Init(new_dim_, num_of_columns_);
+			for(index_t i=0; i<num_of_rows_; i++) {
+			  for(index_t j=0; j<new_dim_; j++) {
+				  w_mat_.set_(i, j, results.get(0, i*new_dim_+j));
+				}
+			}
+			index_t offset_h=num_of_rows_*new_dim_;
+			for(index_t i=0; i<num_of_columns_; i++) {
+			  for(index_t j=0; j<new_dim_; j++) {
+				   h_mat.set(j, i , results(0, offset_h+i*new_dim_+j ));
+				}
+			}
+	 }
+	 void GetW(Matrix *w_mat) {
+	   h_mat->Copy(*w_mat);
+	 }
+	 void GetH(Matrix *h_mat) {
+	   w_mat->Copy(*h_mat);
+	 }
+	 
+ private:
+	fx_module *module_;
+  LBfgs<BigSdpNmfObjective> engine_;
+	BigSdpNmfObjective opt_function_;
+	ArrayList<index_t> rows_;
+	ArrayList<index_t> columns_;
+	ArrayList<double> values_;
+	index_t new_dim_;
+	index_t sdp_rank_;
+	Matrix w_mat_;
+	Matrix h_mat_;
+	index_t num_of_rows_; // number of unique rows, otherwise the size of W
+	index_t num_of_columns_; // number of unique columns, otherwise the size of H
+ 
+	void PreprocessData(Matrix &data_mat) {
+	  values_.Init();
+		rows_.Init();
+		columns_.Init();
+		for(index_t i=0; i<data_mat.n_rows(); i++) {
+		  for(index_t j=0; j< data_mat.n_cols(); j++) {
+			  values_.PushBack(values_.get(i, j));
+				rows_.PushBack(i);
+				columns_.PushBack(j);
+			}
+		}
+	}
+};
+
+#endif
diff --git a/fastlib2/contrib/nvasil/convex_nmf/nmf_objectives.h b/fastlib2/contrib/nvasil/convex_nmf/nmf_objectives.h
index 680853bc65..38c572c94b 100644
--- a/fastlib2/contrib/nvasil/convex_nmf/nmf_objectives.h
+++ b/fastlib2/contrib/nvasil/convex_nmf/nmf_objectives.h
@@ -5,7 +5,7 @@
 
 class BigSdpNmfObjective {
  public:
-	void Init(datanode *module, 
+	void Init(fx_module *module, 
 			ArrayList<index_t> &rows,
 			ArrayList<index_t> &columns,
       ArrayList<double>  &values);
@@ -17,10 +17,11 @@ class BigSdpNmfObjective {
   void UpdateLagrangeMult(Matrix &coordinates);
   void Project(Matrix *coordinates);
   void set_sigma(double sigma); 
-  bool IsDiverging(double objective); 
-
+  void GiveInitMatrix(Matrix *init_data);
+	bool IsDiverging(double objective); 
+ 
  private:
-  datanode *module_;
+  fx_module *module_;
   double sigma_;
 	index_t num_of_columns_;
 	index_t num_of_rows_;
diff --git a/fastlib2/contrib/nvasil/convex_nmf/nmf_objectives_impl.h b/fastlib2/contrib/nvasil/convex_nmf/nmf_objectives_impl.h
index 981ef29278..af224c724d 100644
--- a/fastlib2/contrib/nvasil/convex_nmf/nmf_objectives_impl.h
+++ b/fastlib2/contrib/nvasil/convex_nmf/nmf_objectives_impl.h
@@ -1,5 +1,5 @@
 
-void BigSdpNmfObjective::Init(datanode *module, 
+void BigSdpNmfObjective::Init(fx_module *module, 
 			ArrayList<index_t> &rows,
 			ArrayList<index_t> &columns,
       ArrayList<double>  &values) {
@@ -10,6 +10,8 @@ void BigSdpNmfObjective::Init(datanode *module,
 	num_of_rows_=std::max_element(rows_.begin(), rows_.end())+1;
 	num_of_columns_=std::max_element(columns_.begin(), columns_.end())+1;
 	eq_lagrange_mult_.Init(values_.size());
+  rank_=fx_param_int(module_, "rank", 3);
+	new_dim=fx_param_int(module_, "new_dim", 5);
 }
 
 void BigSdpNmfObjective::Destruct() {
@@ -92,7 +94,18 @@ void BigSdpNmfObjective::Project(Matrix *coordinates) {
 
 void BigSdpNmfObjective::set_sigma(double sigma) {
   sigma_=sigma;
-} 
+}
+
+void BigSdpNmfObjective::GiveInitMatrix(Matrix *init_data) {
+  init_data->Init(rank_, (num_of_rows_+num_of_columns_)*new_dim);
+	for(index_t i=0; i<init_data.n_rows(); i++) {
+	  for(index_t j=0; j<init_data.n_cols(); j++) {
+		  init_data.set(i, j, math::Random());
+		}
+	}
+}
+
+
 bool BigSdpNmfObjective::IsDiverging(double objective) {
   return false;
 }