resistor-cube (empty) → 0.0
raw patch · 4 files changed
+178/−0 lines, 4 filesdep +basedep +hmatrixdep +transformerssetup-changed
Dependencies added: base, hmatrix, transformers, utility-ht
Files
- LICENSE +30/−0
- Setup.hs +2/−0
- resistor-cube.cabal +36/−0
- src/Main.hs +110/−0
+ LICENSE view
@@ -0,0 +1,30 @@+Copyright (c) 2015, Henning Thielemann++All rights reserved.++Redistribution and use in source and binary forms, with or without+modification, are permitted provided that the following conditions are met:++ * Redistributions of source code must retain the above copyright+ notice, this list of conditions and the following disclaimer.++ * Redistributions in binary form must reproduce the above+ copyright notice, this list of conditions and the following+ disclaimer in the documentation and/or other materials provided+ with the distribution.++ * Neither the name of Henning Thielemann nor the names of other+ contributors may be used to endorse or promote products derived+ from this software without specific prior written permission.++THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS+"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT+LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR+A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT+OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,+DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY+THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT+(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ Setup.hs view
@@ -0,0 +1,2 @@+import Distribution.Simple+main = defaultMain
+ resistor-cube.cabal view
@@ -0,0 +1,36 @@+Name: resistor-cube+Version: 0.0+Synopsis: Compute total resistance of a cube of resistors+Description:+ This is an example of how to compute the total resistance+ of a non-trivial circuit of resistors.+ It demonstrates how to build the necessary matrix.+ The computed voltages and currents+ are elements of the null vector of that matrix.+Homepage: http://hub.darcs.net/thielema/resistor-cube+License: BSD3+License-File: LICENSE+Author: Henning Thielemann+Maintainer: haskell@henning-thielemann.de+Category: Math+Build-Type: Simple+Cabal-Version: >=1.10++Source-Repository this+ Tag: 0.0+ Type: darcs+ Location: http://hub.darcs.net/thielema/resistor-cube++Source-Repository head+ Type: darcs+ Location: http://hub.darcs.net/thielema/resistor-cube++Executable resistor-cube+ Main-is: Main.hs+ Build-Depends:+ hmatrix >=0.16 && <0.17,+ transformers >=0.3 && <0.4,+ utility-ht >=0.0 && <0.1,+ base >=4.5 && <4.6+ Hs-Source-Dirs: src+ Default-Language: Haskell2010
+ src/Main.hs view
@@ -0,0 +1,110 @@+{- |+Consider a cube of resistors of equal resistance.+What is the overall resistance from one corner to the opposite one?+-}+module Main where++import qualified Numeric.Container as NC+import qualified Data.Packed.Matrix as Matrix+import qualified Data.Packed.Vector as Vector+import qualified Numeric.LinearAlgebra.HMatrix as HMatrix+import Data.Packed.Matrix (Matrix)++import Control.Applicative (liftA3)+import Control.Functor.HT (outerProduct)+++{-+Set resistance of a primitive resistor to 1.+This way, voltage equals current.+-}+data Coord = C0 | C1 deriving (Eq, Ord, Show, Enum, Bounded)+data Dim = D0 | D1 | D2 deriving (Eq, Ord, Show, Enum, Bounded)+data Corner = Corner Coord Coord Coord deriving (Eq, Ord, Show)+data Edge = Edge Dim Coord Coord deriving (Eq, Ord, Show)++allCoords :: [Coord]+allCoords = [minBound .. maxBound]++allDims :: [Dim]+allDims = [minBound .. maxBound]++flattenCornerIndex :: Corner -> Int+flattenCornerIndex (Corner x y z) =+ (fromEnum x * 2 + fromEnum y) * 2 + fromEnum z++flattenEdgeIndex :: Edge -> Int+flattenEdgeIndex (Edge d x y) =+ (fromEnum d * 2 + fromEnum x) * 2 + fromEnum y+++voltageMatrix :: Matrix Double+voltageMatrix =+ Matrix.fromLists $+ outerProduct+ (\(Edge ed ex ey) c ->+ let ((cx, cy), cz) = selectCornerCoords ed c+ in if ex==cx && ey==cy+ then+ case cz of+ C0 -> 1+ C1 -> -1+ else 0)+ (liftA3 Edge allDims allCoords allCoords)+ (liftA3 Corner allCoords allCoords allCoords)+++-- ToDo: How about cyclic arrangement of dimensions?+selectCornerCoords :: Dim -> Corner -> ((Coord, Coord), Coord)+selectCornerCoords ed (Corner cx cy cz) =+ case ed of+ D0 -> ((cy, cz), cx)+ D1 -> ((cx, cz), cy)+ D2 -> ((cx, cy), cz)++sourceCorner, destCorner :: Corner+sourceCorner = Corner C0 C0 C0+destCorner = Corner C1 C1 C1++currentMatrix :: Matrix Double+currentMatrix =+ Matrix.fromLists $+ outerProduct+ (\c@(Corner _ _ _) (Edge ed ex ey) ->+ let ((cx, cy), cz) = selectCornerCoords ed c+ in if ex==cx && ey==cy+ then+ case cz of+ C0 -> 1+ C1 -> -1+ else 0)+ (filter (/= sourceCorner) $+ filter (/= destCorner) $+ liftA3 Corner allCoords allCoords allCoords)+ (liftA3 Edge allDims allCoords allCoords)+++fullMatrix :: Matrix Double+fullMatrix =+ Matrix.fromBlocks+ [[NC.konst 0 (1,12), Matrix.asRow $ Vector.fromList $ 1 : replicate 7 0],+ [HMatrix.ident 12, voltageMatrix],+ [currentMatrix, NC.konst 0 (6,8)]]++main :: IO ()+main = do+ print fullMatrix+ let [currentVec, potentialVec] =+ Vector.takesV [12,8] $ HMatrix.null1 fullMatrix+ let totalCurrent =+ sum $+ map+ (\d -> NC.atIndex currentVec $ flattenEdgeIndex $ Edge d C0 C0)+ [minBound .. maxBound]+ let cornerPot c = NC.atIndex potentialVec (flattenCornerIndex c)+ let totalVoltage = cornerPot destCorner - cornerPot sourceCorner+ print $ totalVoltage / totalCurrent++{-+result: total resistance is 5/(2+2+2)+-}