Fix build without NPSOL
[openmx:openmx.git] / R / MxCompute.R
1 #
2 #   Copyright 2013 The OpenMx Project
3 #
4 #   Licensed under the Apache License, Version 2.0 (the "License");
5 #   you may not use this file except in compliance with the License.
6 #   You may obtain a copy of the License at
7
8 #        http://www.apache.org/licenses/LICENSE-2.0
9
10 #   Unless required by applicable law or agreed to in writing, software
11 #   distributed under the License is distributed on an "AS IS" BASIS,
12 #   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 #   See the License for the specific language governing permissions and
14 #   limitations under the License.
15
16 setClass(Class = "MxBaseCompute", 
17          representation = representation(
18            id = "integer",
19            "VIRTUAL"),
20          contains = "MxBaseNamed")
21
22 setClassUnion("MxCompute", c("NULL", "MxBaseCompute"))
23
24 setGeneric("convertForBackend",
25         function(.Object, flatModel, model) {
26                 return(standardGeneric("convertForBackend"))
27         })
28
29 setGeneric("assignId",
30         function(.Object, id) {
31                 return(standardGeneric("assignId"))
32         })
33
34 setMethod("assignId", signature("MxBaseCompute"),
35         function(.Object, id) {
36                 .Object@id <- id
37                 .Object
38         })
39
40 setGeneric("getFreeVarGroup",
41         function(.Object) {
42                 return(standardGeneric("getFreeVarGroup"))
43         })
44
45 setMethod("getFreeVarGroup", signature("MxBaseCompute"),
46         function(.Object) {
47                 list()
48         })
49
50 #----------------------------------------------------
51
52 setClass(Class = "MxComputeOperation",
53          contains = "MxBaseCompute",
54          representation = representation(
55            free.set = "MxOptionalChar"))
56
57 setMethod("qualifyNames", signature("MxComputeOperation"),
58         function(.Object, modelname, namespace) {
59                 .Object@name <- imxIdentifier(modelname, .Object@name)
60                 .Object@free.set <- imxConvertIdentifier(.Object@free.set, modelname, namespace)
61                 .Object
62         })
63
64 setMethod("getFreeVarGroup", signature("MxComputeOperation"),
65         function(.Object) {
66                 if (length(.Object@free.set)) {
67                         list(.Object@id, .Object@free.set)
68                 } else {
69                         list()
70                 }
71         })
72
73 setMethod("convertForBackend", signature("MxComputeOperation"),
74         function(.Object, flatModel, model) {
75                 name <- .Object@name
76                 .Object
77         })
78
79 #----------------------------------------------------
80
81 setClass(Class = "MxComputeOnce",
82          contains = "MxComputeOperation",
83          representation = representation(
84            what = "MxCharOrNumber",
85            verbose = "integer",
86            context = "character",
87            maxAbsChange = "logical",
88            fit = "logical",
89            gradient = "logical",
90            hessian = "logical",
91            ihessian = "logical",
92            hgprod = "logical"))
93
94 setMethod("qualifyNames", signature("MxComputeOnce"),
95         function(.Object, modelname, namespace) {
96                 .Object <- callNextMethod();
97                 for (sl in c('what')) {
98                         slot(.Object, sl) <- imxConvertIdentifier(slot(.Object, sl), modelname, namespace)
99                 }
100                 .Object
101         })
102
103 setMethod("convertForBackend", signature("MxComputeOnce"),
104         function(.Object, flatModel, model) {
105                 .Object <- callNextMethod();
106                 name <- .Object@name
107                 if (any(!is.integer(.Object@what))) {
108                         expNum <- match(.Object@what, names(flatModel@expectations))
109                         algNum <- match(.Object@what, append(names(flatModel@algebras),
110                                                              names(flatModel@fitfunctions)))
111                         if (any(is.na(expNum)) && any(is.na(algNum))) {
112                                 stop(paste("Can only apply MxComputeOnce to MxAlgebra or MxExpectation not",
113                                            deparse(.Object@what)))
114                         }
115                         if (!any(is.na(expNum))) {
116                                         # Usually negative numbers indicate matrices; not here
117                                 .Object@what <- - expNum
118                         } else {
119                                 if (any(algNum > length(flatModel@algebras)) && length(algNum) > 1) {
120                                         stop("MxComputeOnce cannot evaluate more than 1 fit function")
121                                 }
122                                 .Object@what <- algNum - 1L
123                         }
124                 }
125                 if (length(.Object@what) == 0) warning("MxComputeOnce with nothing will have no effect")
126                 if (all(.Object@what >= 0) && !.Object@maxAbsChange && !.Object@fit && !.Object@gradient &&
127                             !.Object@hessian && !.Object@ihessian && !.Object@hgprod) {
128                         warning("MxComputeOnce with no action")
129                 }
130                 .Object
131         })
132
133 setMethod("initialize", "MxComputeOnce",
134           function(.Object, what, free.set, context, maxAbsChange, fit, gradient,
135                    hessian, ihessian, hgprod, verbose) {
136                   .Object@name <- 'compute'
137                   .Object@what <- what
138                   .Object@verbose = verbose
139                   .Object@free.set <- free.set
140                   .Object@context <- context
141                   .Object@maxAbsChange <- maxAbsChange
142                   .Object@fit <- fit
143                   .Object@gradient <- gradient
144                   .Object@hessian <- hessian
145                   .Object@ihessian <- ihessian
146                   .Object@hgprod <- hgprod
147                   .Object
148           })
149
150 mxComputeOnce <- function(what, free.set=NULL, context=character(0),
151                           maxAbsChange=FALSE, fit=FALSE, gradient=FALSE,
152                           hessian=FALSE, ihessian=FALSE, hgprod=FALSE, verbose=0L) {
153         new("MxComputeOnce", what, free.set, context, maxAbsChange, fit, gradient,
154             hessian, ihessian, hgprod, verbose)
155 }
156
157 #----------------------------------------------------
158
159 setClass(Class = "MxComputeGradientDescent",
160          contains = "MxComputeOperation",
161          representation = representation(
162            useGradient = "MxOptionalLogical",
163            fitfunction = "MxCharOrNumber",
164            engine = "character",
165            verbose = "integer"))
166
167 setMethod("qualifyNames", signature("MxComputeGradientDescent"),
168         function(.Object, modelname, namespace) {
169                 .Object <- callNextMethod();
170                 for (sl in c('fitfunction')) {
171                         slot(.Object, sl) <- imxConvertIdentifier(slot(.Object, sl), modelname, namespace)
172                 }
173                 .Object
174         })
175
176 setMethod("convertForBackend", signature("MxComputeGradientDescent"),
177         function(.Object, flatModel, model) {
178                 .Object <- callNextMethod();
179                 name <- .Object@name
180                 if (is.character(.Object@fitfunction)) {
181                         .Object@fitfunction <- imxLocateIndex(flatModel, .Object@fitfunction, name)
182                 }
183                 .Object
184         })
185
186 setMethod("initialize", "MxComputeGradientDescent",
187           function(.Object, free.set, engine, fit, useGradient, verbose) {
188                   .Object@name <- 'compute'
189                   .Object@free.set <- free.set
190                   .Object@fitfunction <- fit
191                   .Object@engine <- engine
192                   .Object@useGradient <- useGradient
193                   .Object@verbose <- verbose
194                   .Object
195           })
196
197 mxComputeGradientDescent <- function(type=NULL, free.set=NULL, useGradient=NULL,
198                                      engine=NULL, fitfunction='fitfunction', verbose=0L) {
199 # What to do with 'type'?
200 #       if (length(type) != 1) stop("Specific 1 compute type")
201
202         if (missing(engine)) {
203                 engine <- Sys.getenv("IMX_OPT_ENGINE")
204                 if (!nchar(engine)) {
205                         if (.Call(imxHasNPSOL)) {
206                                 engine <- "NPSOL"
207                         } else {
208                                 engine <- "CSOLNP"
209                         }
210                 }
211         }
212
213         new("MxComputeGradientDescent", free.set, engine, fitfunction, useGradient, verbose)
214 }
215
216 #----------------------------------------------------
217
218 setClass(Class = "MxComputeNewtonRaphson",
219          contains = "MxComputeOperation",
220          representation = representation(
221            fitfunction = "MxCharOrNumber",
222            maxIter = "integer",
223            tolerance = "numeric",
224            verbose = "integer",
225            carefully = "logical"))
226
227 setMethod("qualifyNames", signature("MxComputeNewtonRaphson"),
228         function(.Object, modelname, namespace) {
229                 .Object <- callNextMethod();
230                 for (sl in c('fitfunction')) {
231                         slot(.Object, sl) <- imxConvertIdentifier(slot(.Object, sl), modelname, namespace)
232                 }
233                 .Object
234         })
235
236 setMethod("convertForBackend", signature("MxComputeNewtonRaphson"),
237         function(.Object, flatModel, model) {
238                 .Object <- callNextMethod();
239                 name <- .Object@name
240                 if (is.character(.Object@fitfunction)) {
241                         .Object@fitfunction <- imxLocateIndex(flatModel, .Object@fitfunction, name)
242                 }
243                 .Object
244         })
245
246 setMethod("initialize", "MxComputeNewtonRaphson",
247           function(.Object, free.set, fit, maxIter, tolerance, verbose, carefully) {
248                   .Object@name <- 'compute'
249                   .Object@free.set <- free.set
250                   .Object@fitfunction <- fit
251                   .Object@maxIter <- maxIter
252                   .Object@tolerance <- tolerance
253                   .Object@verbose <- verbose
254                   .Object@carefully <- carefully
255                   .Object
256           })
257
258 mxComputeNewtonRaphson <- function(type, free.set=NULL,
259                                    fitfunction='fitfunction', maxIter = 100L, tolerance=1e-7,
260                                    verbose=0L, carefully=FALSE) {
261
262         new("MxComputeNewtonRaphson", free.set, fitfunction, maxIter, tolerance, verbose, carefully)
263 }
264
265 #----------------------------------------------------
266
267 setClass(Class = "MxComputeSteps",
268          contains = "MxBaseCompute",
269          representation = representation(
270            steps = "list"))
271
272 setMethod("getFreeVarGroup", signature("MxComputeSteps"),
273         function(.Object) {
274                 result <- list()
275                 for (step in .Object@steps) {
276                         got <- getFreeVarGroup(step)
277                         if (length(got)) result <- append(result, got)
278                 }
279                 result
280         })
281
282 setMethod("assignId", signature("MxComputeSteps"),
283         function(.Object, id) {
284                 steps <- .Object@steps
285                 for (sx in 1:length(steps)) {
286                         steps[[sx]] <- assignId(steps[[sx]], id)
287                         id <- steps[[sx]]@id + 1L
288                 }
289                 .Object@steps <- steps
290                 .Object@id <- id
291                 .Object
292         })
293
294 setMethod("qualifyNames", signature("MxComputeSteps"),
295         function(.Object, modelname, namespace) {
296                 .Object@name <- imxIdentifier(modelname, .Object@name)
297                 .Object@steps <- lapply(.Object@steps, function (c) qualifyNames(c, modelname, namespace))
298                 .Object
299         })
300
301 setMethod("convertForBackend", signature("MxComputeSteps"),
302         function(.Object, flatModel, model) {
303                 .Object@steps <- lapply(.Object@steps, function (c) convertForBackend(c, flatModel, model))
304                 .Object
305         })
306
307 #----------------------------------------------------
308
309 setClass(Class = "MxComputeIterate",
310          contains = "MxComputeSteps",
311          representation = representation(
312            maxIter = "integer",
313            tolerance = "numeric",
314            verbose = "integer"))
315
316 setMethod("initialize", "MxComputeIterate",
317           function(.Object, steps, maxIter, tolerance, verbose) {
318                   .Object@name <- 'compute'
319                   .Object@steps <- steps
320                   .Object@maxIter <- maxIter
321                   .Object@tolerance <- tolerance
322                   .Object@verbose <- verbose
323                   .Object
324           })
325
326 mxComputeIterate <- function(steps, maxIter=500L, tolerance=1e-4, verbose=0L) {
327         new("MxComputeIterate", steps=steps, maxIter=maxIter, tolerance=tolerance, verbose)
328 }
329
330 displayMxComputeIterate <- function(opt) {
331         cat(class(opt), omxQuotes(opt@name), '\n')
332         cat("@tolerance :", omxQuotes(opt@tolerance), '\n')
333         cat("@maxIter :", omxQuotes(opt@maxIter), '\n')
334         for (step in 1:length(opt@steps)) {
335                 cat("[[", step, "]] :", class(opt@steps[[step]]), '\n')
336         }
337         invisible(opt)
338 }
339
340 setMethod("print", "MxComputeIterate", function(x, ...) displayMxComputeIterate(x))
341 setMethod("show",  "MxComputeIterate", function(object) displayMxComputeIterate(object))
342
343 #----------------------------------------------------
344
345 setClass(Class = "MxComputeEstimatedHessian",
346          contains = "MxComputeOperation",
347          representation = representation(
348            fitfunction = "MxCharOrNumber",
349            se = "logical"))
350
351 setMethod("qualifyNames", signature("MxComputeEstimatedHessian"),
352         function(.Object, modelname, namespace) {
353                 .Object <- callNextMethod();
354                 .Object@fitfunction <- imxConvertIdentifier(.Object@fitfunction, modelname, namespace)
355                 .Object
356         })
357
358 setMethod("convertForBackend", signature("MxComputeEstimatedHessian"),
359         function(.Object, flatModel, model) {
360                 .Object <- callNextMethod();
361                 name <- .Object@name
362                 if (is.character(.Object@fitfunction)) {
363                         .Object@fitfunction <- imxLocateIndex(flatModel, .Object@fitfunction, name)
364                 }
365                 .Object
366         })
367
368 setMethod("initialize", "MxComputeEstimatedHessian",
369           function(.Object, free.set, fit, want.se) {
370                   .Object@name <- 'compute'
371                   .Object@free.set <- free.set
372                   .Object@fitfunction <- fit
373                   .Object@se <- want.se
374                   .Object
375           })
376
377 mxComputeEstimatedHessian <- function(free.set=NULL, fitfunction='fitfunction', want.se=TRUE) {
378         new("MxComputeEstimatedHessian", free.set, fitfunction, want.se)
379 }
380
381 #----------------------------------------------------
382
383 setClass(Class = "MxComputeSequence",
384          contains = "MxComputeSteps")
385
386 setMethod("initialize", "MxComputeSequence",
387           function(.Object, steps) {
388                   .Object@name <- 'compute'
389                   .Object@steps <- steps
390                   .Object
391           })
392
393 mxComputeSequence <- function(steps) {
394         new("MxComputeSequence", steps=steps)
395 }
396
397 displayMxComputeSequence <- function(opt) {
398         cat(class(opt), omxQuotes(opt@name), '\n')
399         for (step in 1:length(opt@steps)) {
400                 cat("[[", step, "]] :", class(opt@steps[[step]]), '\n')
401         }
402         invisible(opt)
403 }
404
405 setMethod("print", "MxComputeSequence", function(x, ...) displayMxComputeSequence(x))
406 setMethod("show",  "MxComputeSequence", function(object) displayMxComputeSequence(object))
407
408 #----------------------------------------------------
409
410 displayMxComputeOperation <- function(opt) {
411         cat(class(opt), omxQuotes(opt@name), '\n')
412         cat("@id :", opt@id, '\n')
413         cat("@free.set :", omxQuotes(opt@free.set), '\n')
414         invisible(opt)
415 }
416
417 setMethod("print", "MxComputeOperation", function(x, ...) displayMxComputeOperation(x))
418 setMethod("show",  "MxComputeOperation", function(object) displayMxComputeOperation(object))
419
420 displayMxComputeGradientDescent <- function(opt) {
421         cat("@type :", omxQuotes(opt@type), '\n')
422         cat("@engine :", omxQuotes(opt@engine), '\n')
423         cat("@fitfunction :", omxQuotes(opt@fitfunction), '\n')
424         invisible(opt)
425 }
426
427 setMethod("print", "MxComputeGradientDescent",
428           function(x, ...) { callNextMethod(); displayMxComputeGradientDescent(x) })
429 setMethod("show",  "MxComputeGradientDescent",
430           function(object) { callNextMethod(); displayMxComputeGradientDescent(object) })
431
432 convertComputes <- function(flatModel, model) {
433         retval <- lapply(flatModel@computes, function(opt) {
434                 convertForBackend(opt, flatModel, model)
435         })
436         retval
437 }