Upload New File

20210525_Formulation_tool_0.1steps_corrected_int_1_PF_zeroboundary.py

+# -*- coding: utf-8 -*-
+"""
+Product formulator 
+"""
+#Initially all combinations of components are generated. Then all these values are converted in a csv, which can be used in the ternary plot file
+#in R. Later the combinations are subjected to linear regression to obtain all minimal quantity of fractions that adhere to the compositions generated in step 1
+#These values are also converted in a list, which can be copied into the output file that calculates the GWP.
+#In this file only the fractions are used for the formulation window, so without the pea flour. 
+
+
+#Creating all possisbilities for a viscosity range
+import constraint
+import math
+import pandas as pd 
+
+problem = constraint.Problem()
+problem.addVariable('p', [0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1,
+                          1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2,
+                          2.1,2.2,2.3,2.4,2.5,2.6,2.7,2.8,2.9,3,
+                          3.1,3.2,3.3,3.4,3.5,3.6,3.7,3.8,3.9,4, 
+                          4.1,4.2,4.3,4.4,4.5,4.6,4.7,4.8,4.9,5,
+                          5.1,5.2,5.3,5.4,5.5,5.6,5.7,5.8,5.9,6,
+                          6.1,6.2,6.3,6.4,6.5,6.6,6.7,6.8,6.9,7, 
+                          7.1,7.2,7.3,7.4,7.5,7.6,7.7,7.8,7.9,8, 
+                          8.1,8.2,8.3,8.4,8.5,8.6,8.7,8.8,8.9,9, 
+                          9.1,9.2,9.3,9.4,9.5,9.6,9.7,9.8,9.9, 
+                          10.0])
+problem.addVariable('s', [0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1,
+                          1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2,
+                          2.1,2.2,2.3,2.4,2.5,2.6,2.7,2.8,2.9,3,
+                          3.1,3.2,3.3,3.4,3.5,3.6,3.7,3.8,3.9,4, 
+                          4.1,4.2,4.3,4.4,4.5,4.6,4.7,4.8,4.9,5,
+                          5.1,5.2,5.3,5.4,5.5,5.6,5.7,5.8,5.9,6,
+                          6.1,6.2,6.3,6.4,6.5,6.6,6.7,6.8,6.9,7, 
+                          7.1,7.2,7.3,7.4,7.5,7.6,7.7,7.8,7.9,8, 
+                          8.1,8.2,8.3,8.4,8.5,8.6,8.7,8.8,8.9,9, 
+                          9.1,9.2,9.3,9.4,9.5,9.6,9.7,9.8,9.9, 
+                          10.0])
+problem.addVariable('f', [0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1,
+                          1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2,
+                          2.1,2.2,2.3,2.4,2.5,2.6,2.7,2.8,2.9,3,
+                          3.1,3.2,3.3,3.4,3.5,3.6,3.7,3.8,3.9,4,
+                          4.1,4.2,4.3,4.4,4.5,4.6,4.7,4.8,4.9,5,
+                          5.1,5.2,5.3,5.4,5.5,5.6,5.7,5.8,5.9,
+                          6])
+problem.addVariable('r', [0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1,
+                          1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2,
+                          2.1,2.2,2.3,2.4,2.5,2.6,2.7,2.8,2.9,3,
+                          3.1,3.2,3.3,3.4,3.5,3.6,3.7,3.8,3.9,4,
+                          4.1,4.2,4.3,4.4,4.5,4.6,4.7,4.8,4.9,5,
+                          5.1,5.2,5.3,5.4,5.5,5.6,5.7,5.8,5.9,
+                          6])
+
+#Setting desired viscosity add in a range of + and - 2.5%
+var_mean = 1500
+var_boundary = 0.025
+up_vis  = math.log(var_mean * (1 + var_boundary))
+low_vis = math.log(var_mean * (1 - var_boundary))
+
+#Inserting the coefficients of the visocisty regression model 
+def lower_constraint(p, s, f, r):
+    if (p * 0.3064)  + (s * 1.10001) + (f * 1.63015) + (r * 0.45714) + (p * s * -0.03559) + (p * f * -0.0691) + (s * f * -0.22083) >= low_vis:
+        return True
+    
+def upper_constraint(p, s, f, r):
+    if (p * 0.3064)  + (s * 1.10001) + (f * 1.63015) + (r * 0.45714) + (p * s * -0.03559) + (p * f * -0.0691) + (s * f * -0.22083) <= up_vis:
+        return True
+    
+#defining the upper dry matter limit of the training set of the regression model
+def sum_constraint(p, s, f, r):
+    if p + s + f + r <=22:
+        return True
+
+problem.addConstraint(upper_constraint, ['p','s','f','r'])
+problem.addConstraint(lower_constraint, ['p','s','f','r'])
+problem.addConstraint(sum_constraint, ['p','s','f','r'])
+#Obtain dictionnary with all solutions that are loaded into the next step.
+solutions = problem.getSolutions() 
+
+## Making a dataframe of the compositions possible and convert to csv for ternary plot
+SolutionsList_order = [(sub['p'], sub['s'], sub['f'], sub['r']) for sub in solutions] 
+SolutionList = [list(t) for t in SolutionsList_order]
+SolutionList.insert(0,['Protein (DM%)', 'Starch (DM%)', 'Fibre (DM%)', 'Rest (DM%)'])
+df_comp = pd.DataFrame(SolutionList[1:],columns=SolutionList[0])
+df = pd.DataFrame(df_comp.sum(axis = 1))
+df_comp.to_csv (r'C:\ModelsLocal\Functionality models\Ternary plots\1500_0.025mPas_0.1step_22%dm2505_zeroboundary.csv', index = False, header=True)
+
+
+#Creating a 3D scatterplot of all possible combinations of components
+import pandas as pd 
+import plotly.express as px
+from plotly.offline import plot
+
+df_comp = pd.DataFrame(SolutionList[1:],columns=SolutionList[0])
+fig_comp = px.scatter_3d(df_comp, x='Protein (DM%)', y='Starch (DM%)', z='Fibre (DM%)',
+                    color='Rest (DM%)')
+fig_comp.update_layout(
+    font_family="Arial",
+    font_color="black",
+    font_size = 15,
+    )
+
+plot(fig_comp) #online)
+fig_comp.show(renderer = 'png') #in python
+
+
+
+#Convert the found composition into possible combinations of fractions using linear optimisation
+import pulp
+OptimalSolutions = [] #empty list to store the outcome of the variables 
+Fractions = ['CF','FF','FI', 'PF', 'PI', 'SI']
+#Compositions fractions with Eurofins report isolates and Nienkes results MRF
+protein = {'CF': 0.153,
+           'FF': 0.407,
+           'FI': 0.078,
+           'PF': 0.216,
+           'PI': 0.741,
+           'SI': 0.005
+           }
+
+starch = {'CF': 0.626,
+          'FF': 0.138,
+          'FI': 0.338,
+          'PF': 0.479,
+          'PI': 0.004,
+          'SI': 0.914          
+          }
+
+fibre = {'CF': 0.090,
+         'FF': 0.194,
+         'FI': 0.536,
+         'PF': 0.106,
+         'PI': 0.024,
+         'SI': 0.008         
+         }
+
+rest = {'CF': 0.132,
+        'FF': 0.260,
+        'FI': 0.048,
+        'PF': 0.199,
+        'PI': 0.233,
+        'SI': 0.073
+        }
+
+z = 1 #making sure every for loop has a different name 
+for dictionnary in solutions:
+    prob = pulp.LpProblem("The blending problem" + str(z), pulp.LpMinimize)
+    z += 1
+    fractions_vars = pulp.LpVariable.dicts("Frc", Fractions, 0)
+ 
+    #Objective function   
+    prob += pulp.lpSum([fractions_vars[f] for f in Fractions]), "Total fractions in blend" 
+    
+    #Adding constraints
+    prob += pulp.lpSum([fractions_vars[f] for f in Fractions]) <= 22, "Maximum dry matter"
+
+    prob += pulp.lpSum([starch[f] * fractions_vars[f] for f in Fractions]) == (dictionnary["s"]), "StarchRequirement"
+    prob += pulp.lpSum([fibre[f] * fractions_vars[f] for f in Fractions]) == (dictionnary["f"]), "FibreRequirement"
+    prob += pulp.lpSum([rest[f] * fractions_vars[f] for f in Fractions]) == (dictionnary["r"]), "RestRequirement"
+    prob += pulp.lpSum([protein[f] * fractions_vars[f] for f in Fractions]) == (dictionnary["p"]), "ProteinRequirement"
+   
+    #Solving the problem
+    prob.solve()
+    
+    #Saving all solutions in OptimalSolutions
+    if  pulp.LpStatus[prob.status] == "Optimal": #Optimal, Infeasbible
+        i = 1
+        temp_list = []
+        for v in prob.variables():
+            if i < len(prob.variables()):
+                temp_list.append(v.varValue)
+                i += 1
+            elif i == len(prob.variables()):
+                temp_list.append(v.varValue)
+                OptimalSolutions.append(temp_list)
+                temp_list=[]
+                i == 1
+   
+#Transform OptimalSolutions into DF and make figures 
+OptimalSolutions.insert(0,Fractions)
+df_OptimalSolutions = pd.DataFrame(OptimalSolutions[1:], columns = OptimalSolutions[0])
+
+#Create bar chart with all solutions
+ax = df_OptimalSolutions.plot.bar(stacked=True)
+fig_frac = px.bar(df_OptimalSolutions, title="All possible fraction combinations")
+plot(fig_frac) #online)
+fig_frac.show(renderer = 'png') #in python
+
+#Export Optimal Solutions
+df_OptimalSolutions.to_csv (r'C:\Users\lie-p001\.spyder-py3\Product formulator\output_1500_var0.025_averaged_formula_22dm%_PF_CORRECTED_zeroboundary.csv', index = False, header=True)
+
+
+
+
+
+
+
+
+
+
+
+
+