Nothing
#' Derivative calculation of concentration profile
#'
#' Return a the derivative of the concentration profile simulated
#'
#' @param npoints number of points to be used for the derivative
#' @param h space for the finite difference
#' @param Ox data upon the derivative is calculated
#' @param mode "Forward" or "Backward" the derivative will be calculated for the npoints
#' @param Derivative "First" or "Second" derivative to calculate
#' @param CoefMat if T return the derivative coefficient matrix for selected derivative
#'
#'
#' @return a vector with the derivative requested or the coefficient of such derivative
#'
#' @examples
#' Derv(npoints = 2, h = 0.13, Ox = matrix(c(1,2), nrow = 1), mode = "Forward", Derivative = "First")
#'
#' @export
Derv = function(npoints = 2, h, Ox, mode = "Forward", Derivative = "First", CoefMat = FALSE) {
if (mode == "Forward") {
counter = 1
Taylor = ZeroMat(npoints-1, npoints-1)
for (i1 in 1:(npoints-1)){
for (j1 in 1:(npoints-1)) {
Taylor[i1,j1] = (1/factorial(j1))*counter^(j1)
}
counter = counter + 1
}
Ainv = invMat(Taylor)
if (CoefMat == T) {
if (mode == "Forward" & Derivative == "First") {
Coeff = c(-sum(Ainv[1,]), Ainv[1,1:npoints-1])
return(Coeff)
} else if (mode == "Forward" & Derivative == "Second") {
Coeff = c(-sum(Ainv[2,]), Ainv[2,1:npoints-1])
return(Coeff)
}
}
Deriv = matrix(nrow = nrow(Ox), ncol=1)
for (x2 in 1:nrow(Ox)) {
b = c()
for (x1 in 2:npoints){
b = append(b, (Ox[x2,x1] - Ox[x2,1]))
}
if (Derivative == "First") {
Deriv[x2,1] = (1/h)*Ainv[1,]%*%b
} else if (Derivative == "Second") {
Deriv[x2,1] = (1/(h^2))*Ainv[2,]%*%b
}
}
return(Deriv)
} else if (mode == "Backward") {
counter = 1
Taylor = ZeroMat(npoints-1, npoints-1)
for (i1 in 1:(npoints-1)){
for (j1 in 1:(npoints-1)) {
Taylor[i1,j1] = ((-1)^(j1))*(1/factorial(j1))*counter^(j1)
}
counter = counter + 1
}
Ainv = invMat(Taylor)
if (mode == "Backward" & Derivative == "First") {
Coeff = c(-1*Ainv[1,(npoints-1):1], sum(Ainv[1,]))
return(Coeff)
} else if (mode == "Backward" & Derivative == "Second") {
Coeff = c(-1*Ainv[2,(npoints-1):1], sum(Ainv[2,]))
return(Coeff)
}
Deriv = matrix(nrow = nrow(Ox), ncol=1)
for (x2 in 1:nrow(Ox)) {
b = c()
for (x1 in 2:npoints){
b = append(b, (Ox[x2,x1] - Ox[x2,1]))
}
if (Derivative == "First") {
Deriv[x2,1] = (1/h)*Ainv[1,]%*%b
} else if (Derivative == "Second") {
Deriv[x2,1] = (1/(h^2))*Ainv[2,]%*%b
}
}
return(Deriv)
}
}
#' Starting Matrix of oxidazed species
#'
#' Return a matrix ixj filled with 1 value
#'
#' @param i number of rows
#' @param j number of columns
#'
#'
#' @return a matrix of dimention ixj filled with 1 value
#'
#' @examples
#' OneMat(2,2)
#'
#' @export
OneMat = function(i,j=i) {
Y = matrix(data = c(rep(1, i*j)), nrow = i, ncol = j)
return(Y)
}
#' Starting Matrix of reduces species and fluxes
#'
#' Return a matrix ixj filled with 0 value
#'
#' @param i number of rows
#' @param j number of columns
#'
#'
#' @return a matrix of dimention ixj filled with 1 value
#'
#' @examples
#' ZeroMat(2,2)
#'
#' @export
ZeroMat = function(i,j=i){
Y = matrix(data = c(rep(0, i*j)), nrow = i, ncol = j)
return(Y)
}
#' Inverse matrix
#'
#' Returns the inverse matrix of the selected one
#'
#' @param A matrix to be inverted
#'
#'
#' @return inverse matrix of the selected
#'
#' @examples
#' invMat(A = matrix(c(1,2,6,14), nrow = 2))
#'
#' @export
invMat = function(A) {
if (is.matrix(A) == FALSE) {
return("Pointed object is not a matrix")
} else if (det(A) == 0) {
return("Inversion can't be performed on a singular matrix")
}
INV = ZeroMat(nrow(A), ncol(A))
for (i1 in 1:nrow(A)) {
for (j1 in 1:ncol(A)) {
B = A[-i1,-j1]
if (is.matrix(B) == FALSE) {
INV[j1,i1] = (-1)^(i1+j1)*B
} else {
INV[j1,i1] = (-1)^(i1+j1)*det(B)
}
}
}
INV = (1/det(A))*INV
return(INV)
}
#' Parameters call
#'
#' Returns a list with the parameters necessary for the simulation
#'
#' @param Fun Name of the function this function is called to. Must be a string.
#' @param n. Number of electrons
#' @param Temp. Temperature for the simulation
#' @param Dx1. Diffusion coefficient of species One
#' @param eta. OverPotential for potential step
#' @param exptime. experimental time for the simulation
#' @param Eo1. reduction potential of the first electrochemical reaction
#' @param Eo2. reduction potential of the second electrochemical reaction
#' @param Eo3. reduction potential of the third electrochemical reaction
#' @param Eo4. reduction potential of the fourth electrochemical reaction
#' @param Dred1. diffusion coefficient of the first reduced species
#' @param Dred2. diffusion coefficient of the second reduced species
#' @param Dred3. diffusion coefficient of the third reduced species
#' @param Dred4. diffusion coefficient of the fourth reduced species
#' @param alpha1. charge transfer coefficient of the first electrochemical reaction
#' @param alpha2. charge transfer coefficient of the second electrochemical reaction
#' @param alpha3. charge transfer coefficient of the third electrochemical reaction
#' @param alpha4. charge transfer coefficient of the fourth electrochemical reaction
#' @param kc. Chemical rate constant for first Ox Species, used in simulation with just one species
#' @param kco. Chemical rate constant for first Ox Species
#' @param kc1. Chemical rate constant for first Red Species
#' @param ko1. heterogeneous electron transfer rate constant of the first electrochemical reaction
#' @param ko2. heterogeneous electron transfer rate constant of the second electrochemical reaction
#' @param ko3. heterogeneous electron transfer rate constant of the third electrochemical reaction
#' @param ko4. heterogeneous electron transfer rate constant of the fourth electrochemical reaction
#' @param kc2. Chemical rate constant for second Red Species
#' @param kc3. Chemical rate constant for third Red Species
#' @param kc4. Chemical rate constant for fourth Red Species
#' @param Dm. Simulation parameter, maximum 0.5 for explicit methods
#' @param Vi. Initial potential of the sweep
#' @param Vf. Final potential of the sweep
#' @param Vs. Scan rate of the simulation
#' @param l. numer of time steps
#'
#'
#' @return inverse matrix of the selected
#'
#' @examples
#' ParCall("ChronAmp", n. = 1, Temp. = 298, Dx1. = 0.0001, exptime. = 1, Dm. = 0.45, l. = 100)
#'
#' @export
ParCall = function(Fun, n., Temp., Dx1.,
eta., exptime., Eo1., ko1., ko2., kc.,
Dm., Vf., Vi., Vs., alpha1., Eo2., Dred1., Dred2.,
alpha2., Dred3., Dred4., ko3., ko4., kco., kc1., kc2.,
kc3., kc4., alpha3., alpha4., Eo3., Eo4., l.){
if (!(Fun %in% c("ChronAmp", "PotStep", "LinSwp", "CV", "CVEC", "CVEE", "Gen_CV" )) ) {
return("Not suitable function was called for parameter calculation")
}
if (Fun == "ChronAmp") {
FA = 96485
R = 8.3145
f = ((FA*n.)/(R*Temp.))
Da = Dx1./Dx1.
l = l.
tau = exptime.
dt = exptime./l
dtn = dt/tau
h = sqrt((Da*dtn)/Dm.)
j = ceiling(6*(l)^0.5)
vt = c(1:l)
t = dt*vt
Par = list(FA,R,f,dtn,Da,l,h,j,t,tau)
names(Par) = c("FA", "R", "f", "dtn", "Da", "l", "h", "j", "t", "tau")
return(Par)
} else if (Fun == "PotStep") {
FA = 96485
R = 8.3145
f = ((FA*n.)/(R*Temp.))
Da = Dx1./Dx1.
l = l.
tau = exptime.
dt = exptime./l
dtn = dt/tau
p = eta.*f
h = sqrt((Da*dtn)/Dm.)
j = ceiling(6*(l)^0.5)
vt = c(1:l)
t = dt*vt
Par = list(FA,R,f,dtn,p,Da,l,h,j,t,tau)
names(Par) = c("FA", "R", "f", "dtn", "p", "Da", "l", "h", "j", "t","tau")
return(Par)
} else if (Fun == "LinSwp") {
FA = 96485
R = 8.3145
f = ((FA*n.)/(R*Temp.))
Da = Dx1./Dx1.
exptime = abs(Vf.-Vi.)/Vs.
l = l.
dt = exptime/l
tau = (1/(f*Vs.))
dtn = dt/tau
j = ceiling(6*(l)^0.5)
h = sqrt((Da*dtn)/Dm.)
vt = c(1:(l))
t = dt*vt
PotentialScan = Vi.-Vs.*t
p = f*(Vi.- Eo1.) - t/tau
KO = ko1.*sqrt(tau/(Dx1./10000))
Kf = KO*exp(-alpha1.*p)
Kb = KO*exp((1-alpha1.)*p)
Par = list(FA,R,f,dtn,Da,l,h,j,t,PotentialScan,KO,Kf,Kb,tau, p)
names(Par) = c("FA", "R", "f", "dtn", "Da", "l", "h",
"j", "t", "PotentialScan", "KO", "Kf", "Kb", "tau", "p")
return(Par)
} else if (Fun == "CV" | Fun == "CVEC") {
FA = 96485
R = 8.3145
f = ((FA*n.)/(R*Temp.))
Da = Dx1./Dx1.
exptime = 2*abs(Vf.-Vi.)/Vs.
l = l.
dt = exptime/l
tau = (1/(f*Vs.))
dtn = dt/tau
j = ceiling(6*(l)^0.5)
h = sqrt((Da*dtn)/Dm.)
vt = c(1:(l))
t = dt*vt
forwardScan = Vi.-Vs.*t[1:((l/2))]
backwardScan = Vf. + Vs.*t[1:((l/2))]
PotentialScan = c(forwardScan, backwardScan)
pf = f*(Vi.- Eo1.) - t[1:((l/2))]/tau
pb = f*(Vf.- Eo1.) + t[1:((l/2))]/tau
p = c(pf,pb)
KO = ko1.*sqrt(tau/(Dx1./10000))
Kf = KO*exp(-alpha1.*p)
Kb = KO*exp((1-alpha1.)*p)
if (Fun == "CVEC") {
KC = kc.*tau
Par = list(FA,R,f,dtn,Da,l,h,j,t,PotentialScan,KO,Kf,Kb,KC,tau,p)
names(Par) = c("FA", "R", "f", "dtn", "Da", "l", "h",
"j", "t", "PotentialScan", "KO", "Kf", "Kb", "KC","tau", "p")
return(Par)
}
Par = list(FA,R,f,dtn,Da,l,h,j,t,PotentialScan,KO,Kf,Kb,tau,p)
names(Par) = c("FA", "R", "f", "dtn", "Da", "l", "h",
"j", "t", "PotentialScan", "KO", "Kf", "Kb","tau","p")
return(Par)
} else if (Fun == "CVEE"){
FA = 96485
R = 8.3145
f = ((FA*n.)/(R*Temp.))
DOx = Dx1./Dx1.
DRED = Dred1./Dx1.
DRED2 = Dred2./Dx1.
exptime = 2*abs(Vf.-Vi.)/Vs.
l = l.
dt = exptime/l
tau = (1/(f*Vs.))
dtn = dt/tau
j = ceiling(6*(l)^0.5)
h = sqrt((DOx*dtn)/Dm.)
vt = c(1:(l))
t = dt*vt
forwardScan = Vi.-Vs.*t[1:((l/2))]
backwardScan = Vf. + Vs.*t[1:((l/2))]
PotentialScan = c(forwardScan, backwardScan)
p1f = f*(Vi.- Eo1.) - t[1:((l/2) +1)]/tau
p1b = f*(Vf.- Eo1.) + t[1:((l/2) +1)]/tau
p2f = f*(Vi.- Eo2.) - t[1:((l/2) +1)]/tau
p2b = f*(Vf.- Eo2.) + t[1:((l/2) +1)]/tau
p1 = c(p1f,p1b)
p2 = c(p2f,p2b)
KO1 = ko1.*sqrt(tau/(Dx1./10000))
KO2 = ko2.*sqrt(tau/(Dx1./10000))
Kf1 = KO1*exp(-alpha1.*p1)
Kf2 = KO2*exp(-alpha2.*p2)
Kb1 = KO1*exp((1-alpha1.)*p1)
Kb2 = KO2*exp((1-alpha2.)*p2)
Par = list(FA,R,f,dtn,DOx,DRED,DRED2,l,h,j,t,PotentialScan,
KO1,KO2,Kf1,Kf2,Kb1,Kb2,tau,p1,p2)
names(Par) = c("FA", "R", "f", "dtn", "DOx", "DRED", "DRED2", "l", "h",
"j", "t", "PotentialScan", "KO1",
"KO2", "Kf1", "Kf2", "Kb1", "Kb2","tau", "p1", "p2")
return(Par)
} else if (Fun == "Gen_CV"){
FA = 96485
R = 8.3145
f = ((FA*n.)/(R*Temp.))
DOx = Dx1./Dx1.
DRED = Dred1./Dx1.
DRED2 = Dred2./Dx1.
DRED3 = Dred3./Dx1.
DRED4 = Dred4./Dx1.
exptime = 2*abs(Vf.-Vi.)/Vs.
l = l.
dt = exptime/l
tau = (1/(f*Vs.))
dtn = dt/tau
j = ceiling(6*(l)^0.5)
h = sqrt((DOx*dtn)/Dm.)
vt = c(1:(l))
t = dt*vt
forwardScan = Vi.-Vs.*t[1:((l/2))]
backwardScan = Vf. + Vs.*t[1:((l/2))]
PotentialScan = c(forwardScan, backwardScan)
p1f = f*(Vi.- Eo1.) - t[1:((l/2))]/tau
p1b = f*(Vf.- Eo1.) + t[1:((l/2))]/tau
p2f = f*(Vi.- Eo2.) - t[1:((l/2))]/tau
p2b = f*(Vf.- Eo2.) + t[1:((l/2))]/tau
p3f = f*(Vi.- Eo3.) - t[1:((l/2))]/tau
p3b = f*(Vf.- Eo3.) + t[1:((l/2))]/tau
p4f = f*(Vi.- Eo4.) - t[1:((l/2))]/tau
p4b = f*(Vf.- Eo4.) + t[1:((l/2))]/tau
p1 = c(p1f,p1b)
p2 = c(p2f,p2b)
p3 = c(p3f,p3b)
p4 = c(p4f,p4b)
KO1 = ko1.*sqrt(tau/(Dx1./10000))
KO2 = ko2.*sqrt(tau/(Dx1./10000))
KO3 = ko3.*sqrt(tau/(Dx1./10000))
KO4 = ko4.*sqrt(tau/(Dx1./10000))
Kf1 = KO1*exp(-alpha1.*p1)
Kf2 = KO2*exp(-alpha2.*p2)
Kf3 = KO3*exp(-alpha3.*p3)
Kf4 = KO4*exp(-alpha4.*p4)
Kb1 = KO1*exp((1-alpha1.)*p1)
Kb2 = KO2*exp((1-alpha2.)*p2)
Kb3 = KO3*exp((1-alpha3.)*p3)
Kb4 = KO4*exp((1-alpha4.)*p4)
KCo = kco.*tau
KC1 = kc1.*tau
KC2 = kc2.*tau
KC3 = kc3.*tau
KC4 = kc4.*tau
Par = list(FA,R,f,dtn,DOx,DRED,DRED2,DRED3,DRED4,l,h,j,t,PotentialScan,
KO1,KO2,KO3,KO4,Kf1,Kf2,Kf3,Kf4,
Kb1,Kb2,Kb3,Kb4,KCo,KC1,KC2,KC3,KC4,tau,p1,p2,p3,p4)
names(Par) = c("FA", "R", "f", "dtn", "DOx", "DRED", "DRED2", "DRED3",
"DRED4", "l", "h",
"j", "t", "PotentialScan", "KO1",
"KO2", "KO3", "KO4", "Kf1", "Kf2",
"Kf3", "Kf4", "Kb1", "Kb2", "Kb3", "Kb4",
"KCo", "KC1", "KC2", "KC3", "KC4", "tau",
"p1", "p2", "p3", "p4")
return(Par)
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.